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ABSTRACT 


The ability of the U.S. Navy’s Coupled Ocean-Atmosphere Mesoscale Prediction 
System (COAMPS™) to accurately forecast the height and structure of the Marine 
Boundary Layer (MBL) in the coastal zone is analyzed and compared to surface and 
aircraft observations from the Dynamics and Evolution of Coastal Stratus (DECS) field 
study conducted along the central coast of California from June, 16 to July 22, 1999. The 
stratus field was found to have significant mesoscale variability within 100 km of the 
coast due to interaction between the mean flow and the coastal terrain. This structure is 
consistent with general hydraulic flow theory and the development of a low level coastal 
jet; however the specific characteristics on any given day were very sensitive to flow 
direction, inversion height, and synoptic conditions. With some modifications, the model 
predicted the general evolution of these events with qualitative fidelity but was slow to 
dissipate the cloud and frequently produced surface fog vice stratus. 

A consistent tendency was found in the model predictions of inversion heights 
200-300 meters too low, weak inversion strengths, high integrated liquid water content, 
and weak buoyancy flux near the cloud top. These observed biases are consistent with 
underestimating the cloud top entrainment velocity and entrainment fluxes in the 
modeled boundary layer. An explicit entrainment parameterization was developed to 
better represent the sub-grid scale processes at cloud top and tested in the single column 
and three-dimensional versions of COAMPS™. The first step in this process is accurate 
determination of the inversion height. It was found that the current method of 
determining boundary layer height in COAMPS™ based on the bulk Richardson number 
frequently misdiagnosed the boundary layer height as occurring in the subcloud layer 
when a weakly stable surface-based inversion was present. Alternative methods based on 
the liquid water content and the liquid potential temperature gradient were tested and 
showed a more consistent diagnosis for the observed conditions in the marine boundary 
layer. The explicit entrainment parameterization was found to generally improve the 
boundary layer height and cloud liquid water content as compared to field observations, 
however the modeled boundary layer still exhibited a low bias and the entrainment 
velocity was higher than is generally expected from field studies for this regime. 
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I. INTRODUCTION 

A. THE IMPORTANCE OF ACCURATE MARINE BOUNDARY LAYER 

FORECASTS TO THE U.S. NAVY 

The Navy operates primarily in the Marine Atmospheric Boundary Layer 
(MABL). Accurate analysis and forecasting of sensible weather elements, especially 
those that impact aviation such as ceiling and visibility, are critical to safe and efficient 
operations. In addition, the large temperature and moisture gradients found in the MABL 
cause rapid and significant changes in refractive conditions that impact radars, 
electromagnetic and electro-optic surveillance and targeting, and communications. To 
better understand this environment, several recent field experiments have focused on the 
MABL cloud cover and dynamics. These include the First ISCCP Regional Experiment 
(FIRE; Cox et al. 1987), Dynamics and Chemistry of Marine Stratocumulus (DYCOMS I 
and II; Lenschow et al. 1988; Stevens et al. 2003a) in the California Bight, CO AMPS™ 
Operational Satellite and Aircraft Test (COSAT; Wetzel et al. 2001) off the Oregon 
coast, Atlantic Stratocumulus Transition Experiment (ASTEX; Albrecht et al. 1995) in 
the northeast Atlantic near the Azores, and Development and Evolution of Coastal 
Stratocumulus (DECS; Kalogiros and Wang 2002) in Monterey Bay. Despite improved 
understanding of the relevant processes, accurate prediction of the MABL diurnal 
evolution is still limited by the scarcity of meteorological measurements over oceanic 
regions to correctly prescribe initial conditions and the lag between new understanding 
from field data and the appropriate realization of this understanding in numerical 
prediction schemes. 


B. GENERAL MODEL OF THE COASTAL MARINE BOUNDARY LAYER 

The marine boundary layer inversion height in the subtropical Eastern Pacific is 
the result of a balance between large-scale subsidence due to the semi-permanent 
subtropical anticyclones and turbulent transport. Turbulent eddies result in entrainment 
of free atmosphere air into the boundary layer, which can lead to boundary layer growth 
and weakening of the stable inversion at the boundary layer top if the free atmosphere 
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above the entrainment zone is weakly stratified. Because surface sensible heat fluxes are 
small in the region of cold ocean boundary currents, surface shear, surface moisture flux, 
in-cloud latent heat release, and cloud-top radiative cooling largely generate the 
turbulence in the mixed layer. Surface buoyancy flux is primarily generated by moisture 
flux, and sensible heat flux is usually small and may even be negative. Additional factors 
in the coastal zone are less well understood and include the role of upwelling and 
mesoscale variability in the sea surface temperatures, the role of shear at and above cloud 
top due to the presence of a low level coastal jet, the quantitative effects on entrainment 
rate due to evaporative cooling as the cloud layer mixes with dry inversion air near cloud 
top, and the impact of sea breeze and coastal jet circulations on cloud evolution, 
microphysics and drizzle. Furthermore, local divergence patterns can be greatly affected 
by flow interaction with the coastal topography leading to complex and diumally varying 
mesoscale structure in the inversion strength and height. The rapid change in surface 
characteristics at the coast also leads to complex coastal circulations and thermal internal 
boundary layers. This research will focus in particular on the role of cloud-top 
entrainment and its representation in mesoscale models, the interaction between cloud- 
top processes and mechanisms particular to the coastal zone such as coastal jet-related 
vertical shear and on horizontal variability in the forcing and resultant boundary layer 
structure in the coastal region. 


C. RECENT STUDIES IN IMPROVING MESOSCALE MODELING USING 

COAMPS™ 

Several researchers supported by the ONR and other agencies in the Department 
of Defense (DoD) are working on improving the numerical representation in CO AMPS™ 
of various mesoscale aspects of the MABL such as the surface layer parameterization 
(Wang 2002, Thompson 2002), the radiation scheme (Kong 2002), the cloud 
microphysics and drizzle parameterizations (Mechem and Kogan 2003), and the 
turbulence parameterization (Whisenhant 2003). Case studies of the Cloud-Topped 
marine Boundary Layer (CTBL) structure (e.g. Haack and Burk 2001; Burk et al. 2003) 
and the numerical representation of CTBL (Wetzel et al. 2001) show reasonable 


2 



qualitative skill but still do not represent the MABL with sufficient fidelity for 
operational prediction. Kong (2002) specifically investigated the predictive skill of 
COAMPS™ for a point verification of stratus and found skill in predicting dissipation 
driven by short wave radiation, but required very high resolution of 50 m in the vertical 
and 2 km in the horizontal, a modified radiation scheme, and an auto-conversion 
threshold of 0.5 g/ kg to force drizzle processes. Furthermore, none of these 
modifications were evaluated for impact to predictive skill in other geographic areas or 
regimes. Nowcasting efforts combining mesoscale numerical prediction models with 
climatology averaged initial structure nontraditional data sources (Leidner et al. 2001) 
and rules based expert systems (Hilliker and Fritsch 1999), have shown promise but 
require local training observational data sets and are not suitable for Naval use because 
the area of interest in Naval operations changes rapidly. A preferable solution that we 
pursue in this study is to improve the model representation of dynamic processes such as 
cloud-top entrainment in a physically consistent and numerically efficient way. 


D. RESEARCH OBJECTIVES 

The overall objective of this study is to better understand the coastal marine 
stratocumulus evolution and its interaction with the complex flow field resulting from 
interaction with the adjacent terrain. Particular emphasis has been placed on the role of 
cloud-top entrainment in this development. This study theorizes that entrainment of 
warm, dry inversion air into the cloud top driven by cloud-top turbulent processes has a 
significant impact on the mesoscale structure and evolution of the CTBL. This process 
occurs at vertical scales of order 10 meters, which is below typical vertical resolution in 
operational forecast models. It is theorized that better representation of the entrainment 
processes and resulting dynamic fluxes will improve operational mesoscale prediction. 
Wetzel et al. (2001), Haack and Burk (2001) and preliminary data from DECS and 
DYCOMS II have indicated that, at resolutions usable for operational prediction, 
COAMPS™ reproduces the general structure and trends in the CTBL but the marine 
boundary layer is too shallow, has too thick a cloud layer, and has liquid water path 
(LWP) values 3-4 times larger than expected. The forecast model in their study produced 
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reasonable area coverage but Liquid Water Path (LWP) was 3-4 times larger than 
observations derived from satellite data. The observed model bias of high LWP could 
result from either lack of precipitation processes or under-represented cloud top 
entrainment of warm dry inversion air while the shallowness of the layer is likely due to 
under-prediction of the turbulent entrainment rate. Mechem and Kogan (2003) 
demonstrated that high vertical (—50 m) and horizontal (2 km) resolution simulations with 
a 2-moment drizzle scheme (Khairoutdinov and Kogan 2000) can produce a CTBL 
qualitatively similar in vertical structure, diurnal variability, horizontal inhomogeneity to 
an observed case and can reproduce drizzle-induced stratocumulus to cumulus transition. 
Preliminary analysis in this study indicates that, when the auto-conversion threshold for 
cloud liquid water conversion to drizzle is adjusted in the lower levels from the default 
value of 1.0 g kg' 1 suitable for convective rain, to a value of 0.5 g kg' 1 for marine stratus 
as recommended by Mechem and Kogan (2003), the model LWC and surface 
accumulated precipitation amounts better match the observations. Other aspects of the 
boundary layer structure such as inversion height and strength are not affected and the 
LWC is generally 20-25% higher than observations rather than the 300-400% found in 
the COSAT experiment. This study examines if an explicit entrainment parameterization 
which better represents the effects of sub-grid scale (SGS) processes at cloud top will 
result in a modeled CTBL structure that more closely resembles the observed CTBL 
structure. 

In support of the overall objective, this research analyzed the observed boundary 
layer structure during the months of June and July 1999 in the DECS dataset as compared 
to model simulations under a variety of synoptic conditions and mesoscale coastal 
processes. Observed and simulated turbulent fluxes and thermodynamic forcing were 
examined with respect to diurnal evolution and horizontal variability of the CTBL in the 
coastal zone. An explicit entrainment parameterization was implemented following 
Grenier and Bretherton (2001) and the resulting fluxes and boundary layer structure were 
compared to observed and modeled fluxes in the current method used in COAMPS™. 

Chapter II summarizes the relevant background and research papers related to this 
study; Chapter III discusses methodology, observational datasets, and general control run 
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model results; Chapter IV discusses methods of identifying the boundary layer height in 
COAMPS™ and its role in parameterizing entrainment; Chapter V discusses the 
development of an explicit entrainment parameterization and single column model 
(SCM) results; Chapter VI discusses three dimensional results, and Chapter VII 
summarizes our findings and provides recommendations for future research. 
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II. BACKGROUND 


A. CLOUD TOPPED BOUNDARY LAYERS 

Modeling of the CTBL in CO AMPS™ is the main focus of this research, thus it is 
important to understand the physical processes that generate and maintain turbulence and 
result in the observed structure and evolution. The marine atmospheric boundary layer 
(MABL) is the lowest (-10%) part of the troposphere that is directly affected by ocean 
surface forcing on a time scale of order 1 hour or less (Stull 1988). The MABL is 
comprised of the interfacial layer, where molecular diffusion dominates (order ~1 cm), 
the surface layer, (order -10% of MABL height) where turbulent fluxes are nearly 
constant with height, and the outer layer where turbulent fluxes are linear with height if 
well mixed. A balance between turbulent mixing and large-scale subsidence determines 
the depth of the MABL. In the case of the eastern subtropical ocean regions, the 
subtropical high produces a permanent subsidence inversion and generally northerly flow 
over increasing SST. This creates a nearly permanent MABL structure of increasing 
boundary layer heights southward and westward from the coast. In summertime, the 
inversion base is frequently above the lifting condensation level resulting in persistent 
large-scale sheets of stratus and stratocumulus which have important effects on the total 
global radiation budget as well as significant climate and economic impacts on coastal 
communities. Burk and Thompson (1996) and Dorman et al. (2000) provide good 
descriptions of the summertime MABL along the U.S. West coast. 

Stratiform low cloud is a common feature through most of the world’s oceans. 
Marine stratocumulus CTBL in particular are often observed along western continental 
shelves at the eastern edge of the summertime subtropical high. As discussed in Klein 
and Hartmann (1993), from a climate perspective, marine stratiform clouds can be 
classed into three basic types as follows. Subtropical CTBL form primarily in summer 
over cold western boundary currents on the east side of subtropical highs and are 
maintained by a balance between large-scale subsidence and convection driven by cloud 
top radiative cooling. Midlatitude CTBL form in winter as cold continental air flows 
over warm western boundary currents. Arctic stratus forms from radiative cooling of 
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moist air entering the Arctic from the subpolar regions. This study will primarily focus 
on subtropical stratus, although the same boundary layer model can be used to describe 
the other types if parameterizations are kept physically based and consistent. As depicted 
schematically in Figure 2.1 from Atkinson and Zhang (1996), the general mixed layer 
structure in the cloud-topped case is primarily a balance between downward adiabatic 
warming from large-scale subsidence, cloud-top radiative cooling, cloud dynamics and 
thermodynamics, and weak surface forcing. 

As discussed in Paluch and Lenschow (1991) based on FIRE data, the 
stratocumulus CTBL has a more complex structure than the typical clear convective 
boundary layer. This is because of additional mechanisms affecting the turbulence 
structure due to thermodynamic effects (primarily water phase changes), longwave and 
shortwave radiative effects associated with the cloud layer, and heat and moisture 
transport by drizzle falling and re-evaporation. The inversion at the cloud top tends to be 
sharper in both the virtual potential temperature (0 V ) and total water mixing ratio (q t ) than 
in clear boundary layers due to enhanced turbulent mixing in the upper part of the 
boundary layer driven mainly by cloud-top radiative cooling. This longwave radiative 
cooling and shortwave daytime heating modify the in-cloud thermodynamic profile and 
produce turbulence, although the net effect is also modified by evaporative cooling or 
sensible heating as warm, dry inversion air is entrained and mixed with cloudy parcels. 
Some researchers have proposed that this entrainment mixing can lead to buoyancy 
reversal and cloud dissipation through Cloud Top Entrainment Instability (CTEI) 
(Randall 1980, Deardorff 1980) if the buoyancy jump condition is sufficiently large, 
although the exact criteria for CTEI to occur is still a subject of active research. In the 
absence of subsidence, entrainment will tend to deepen the boundary layer by mixing free 
atmosphere air into the layer and increasing the mass. In-cloud radiative warming or sub¬ 
cloud evaporative cooling due to drizzle can also lead to less well-mixed or “decoupled” 
boundary layers which can lead to cloud dissipation or pooling of heat and moisture in 
the surface layer and shallow cumulus convection. 

Drizzle is frequently observed in CTBL. As discussed in Jiang et al. (2002), the 
role of drizzle in CTBL structure is quite variable and can lead to cloud base lowering, 
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cloud dissipation, or stratocumulus to cumulus transition. When drizzle reaches the 
surface, the entire sub-cloud layer is evaporatively cooled and moistened, resulting in 
reduced surface buoyancy flux. This can increase the surface temperature gradient and 
result in transition from stratus to shallow convection and cumulus transition. For cases 
where the precipitation remains aloft (virga), the effect is variable, with lightly drizzling 
(-0.1 mm day’ 1 ) cases becoming stabilized and decoupled from the surface, possibly 
leading to daytime dissipation, and heavier rates (~0.5 mm day-1) increasing sub-cloud 
instability thus increasing development of penetrating cumulus, turbulence generation, 
and entrainment (Wang and Wang 1994). These effects are also dependent on the large- 
scale subsidence and the aerosol loading (e.g. Twomey 1977, Albrecht 1989, Penner et al. 
2001). 
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Figure 2.1. Primary physical processes in StCu topped boundary layers. From 
Atkinson and Zhang (1996). 
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B. THE SUMMERTIME CTBL STRUCTURE ALONG THE CALIFORNIA 

COAST 

1. The Low Level Coastal Jet (LLCJ) 

The California coastal region in particular has been studied extensively through 

both field observations and modeling studies as a particular case of cloud-topped 

boundary layers. As described in Beardsley et al. (1987) and Burk and Thompson 

(1996), a primary feature of the summertime California coastal region is a persistent, 

coast-parallel northerly low level jet (LLCJ) that forms due to the cross-coast pressure 

gradient between the East Pacific subtropical high, Sonoran Desert thermal low, capping 

inversion, and coastal mountain range. As depicted in Figures 2.2 and 2.3, the horizontal 

wind speed increases with decreasing height and coastward distance, as would be 

expected from thermal wind considerations, to a maximum at the top of the boundary 

layer inversion, and then decreases due to turbulent diffusion to the surface. This results 

in an inversion that slopes more sharply in the inner coastal zone than the gradual 

westward and southward slope expected from synoptic level considerations and a distinct 

diurnal character due to additional diurnal temperature fluctuations over land. This 

structure is well represented in COAMPS™, as can be seen in Figure 2.4 of day and night 

East-West cross-sections of the boundary layer at 37.1 N, which is across the southern 

part of Monterey Bay. As discussed in Burk and Thompson (1996), Cui et al. (1998), 

Dorman et al. (2000) and elsewhere, the exact mechanisms causing the boundary layer 

and jet structure and spatial and temporal variability are still areas of active research, 

however a proper specification of the coastal mountains and the Sierra Nevada are critical 

to the structure of the MABL and LLCJ close to the coast (Cui et al. 1998). The dynamic 

impact of the marine stratocumulus on the wind field has received little attention, but is 

potentially an important factor as it’s a common characteristic in this regime. In addition 

to the broader jet feature that often extends to 200-300 km from shore, the inner coastal 

region, defined here as less than 100 km off shore, has an additional mesoscale structure 

of relative maximums and minimums in wind speed and boundary layer height 

respectively downstream of prominent points. These are associated with regions of 

increased horizontal divergence and forced subsidence as depicted in Figure 2.5. The 

reverse structure is also seen upstream of these points. The four dimensional structure of 
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these features is difficult to observe directly but is well resolved in operational mesoscale 
models as discussed in Burk and Thompson (1996) and Koracin and Dorman (2001). 
The model simulations show reasonable agreement with surface observations in 
horizontal extent, magnitude, and diurnal variability with the observed surface wind 
maximum close to 0500 UTC, just after sunset. Based on mesoscale model simulations 
and observations, Koracin and Dorman (2001) found the mesoscale low level coastal jet 
(MLLCJ), embedded within the broader LLCJ feature, can largely be described by 
shallow water hydraulic flow theory, with the mesoscale variability driven by flow along 
a complex coastline constrained by the thermally induced cross-coast gradient. Bends in 
the coast result in low wind, high boundary layer compression bulges upstream from 
points, high wind, low boundary layer expansion fans downstream in the capes, and 
transcritical flow further downstream where the flow recovers to the gradient value. 
Although it’s been noticed that the persistence or dissipation of coastal stratus is often 
correlated with this mesoscale divergence and vertical wind shear pattern, there is not 
much in the literature about the exact mechanisms involved. 
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Figure 2.2. Basic structure of the CTBL along the California coast. From 
Beardsley et al. (1987). 
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Figure 2.3. Diurnal variability in the inversion structure and vertical wind profile 
along the summertime California coast. From Beardsley et al. (1987). 
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Figure 2.4. Simulated diurnal virtual potential temperature (K) and liquid water 
mixing ratio (g kg' 1 ) contours in an East-West cross-section from Monterey Bay, 
CA (36.71 N) (between Pigeon Point and Point Sur) for July 6-7, 1999, a day with 
a strong LLCJ. Note the a) daytime shallowing and b) nighttime deepening of the 
boundary layer at approximately 122 W consistent with the convergence/ 
divergence pattern in Figure 2.4. 
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Figure 2.5. Simulated diurnal boundary layer height and divergence pattern along 
the California coast. Panels a) and c) are daytime (1400 LST) and b) and d) are 
nighttime (0200 LST). Adapted from Koracin and Dorman (2001). 
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2. The Sea Breeze/ Mountain Valley Circulation (SBMV) 


In addition to a spatially and diumally varying LLCJ, several points along the 
California coast develop a robust sea breeze mountain valley circulation (SBMV) within 
approximately 100 km of the coast with an increasing diurnal amplitude shoreward 
(Kindle et al. 2002). This circulation is particularly prominent in breaks in the coastal 
mountains such as at San Francisco, Monterey Bay (Banta et al. 1993) and the Los 
Angeles Basin (Lu and Turco 1994). The SBMV differs in timing and structure from a 
classic sea breeze front in that upslope flow and elevated solar driven heating in complex 
terrain can change the timing and orientation of the flow. Lidar studies by Banta et al. 
(1993) of the Monterey Bay showed that the SBMV grew during the late morning hours 
similar to other sea breeze regions, but grew much faster over land than over the cold 
coastal waters. The SBMV also varied from a classic sea breeze in that there was little 
change in vertical velocity at the sea breeze front due to strong synoptic subsidence 
extending to at least 1 km, and there was little evidence of a return flow aloft towards the 
ocean. Horizontal extent seaward seems to be on the order of 30-50 km, which is also 
supported by model and buoy data as discussed in Kindle et al. (2002) using buoy data 
and CO AMPS™ simulations to explore the surface wind short and long period variability 
over Monterey Bay. The subsiding return flow of this sea breeze circulation contributes 
to the observed lowering of the inversion height near the coast in late afternoon but does 
not independently cause the observed pattern as the convergence pattern remains robust 
over a longer diurnal period and is also observed downstream of points such as Point Sur 
which do not develop a clear SBMV. 

3. Coastal Fog and Stratus 

Several studies of fog and stratus have been done in Northern California due to 
the critical operational and economic impact low ceilings and visibility has on marine and 
aviation communities. The impact on flight scheduling at San Francisco International 
airport (SFO) is discussed in Hilliker and Fritsch (1999) which, as seen in Figure 2.6, 
found a relative maximum in occurrence of low stratus and fog at SFO at 1500 UTC (just 

after sunrise) and the greatest decrease in frequency of occurrence 1-4 hours afterwards, 
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indicating the most likely time of dissipation leading to a relative minimum at 2300 UTC 
(mid afternoon). Figure 2.7 depicts a typical diurnal variation in the large-scale 
distribution. From a prediction perspective, few studies have attempted to distinguish 
between the dominant physical mechanisms underlying the evolution of the inner coastal 
stratus, inner coastal surface fog, outer coastal stratus, and valley fog. Koracin et al. 
(2001) examined the transition between stratus and fog over water from a Lagrangian 
perspective. They determined that fog forms in response to a relatively long 
preconditioning of the marine boundary layer, with cloud top cooling as the primary 
mixing mechanism for the lifting of the inversion resulting in fog developing from 
stratus, and subsidence as the predominant cause of lowering of the inversion base 
resulting in stratus transitioning to fog. As will be discussed in Chapter III, it was found 
that in order to properly model the CTBL, the long preconditioning required for the 
CTBL to develop and reach semi-equilibrium requires an upstream model domain on the 
order of 200-300 km and a 36-48 hour data assimilation period (“warm start”) to 
accurately depict the CTBL. 
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Figure 2.6. Diurnal occurrence of a) ceilings less than 914 m at SFO from May to 
October and b) corresponding hourly arrival demand of air traffic. From Hilliker 
and Fritsch (1999). 
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Figure 2.7. Typical diurnal pattern in summertime coastal cloud distribution. 
Note nocturnal development and daytime dissipation of valley fog (VF), terrain 
forced coastal clearing/ persistence around points (TFV) and broader, synoptically- 
driven clearing (SynC) over water that still has a coastal characteristic. Figures a) 
and c) are early morning at 1530 UTC and figure b) is at 0130 UTC late afternoon. 
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C. REPRESENTATION OF BOUNDARY LAYER TURBULENCE IN 
COAMPS™ 


In numerical models, the choice of parameterization scheme for turbulent flux 
calculation is restricted by available computer power and time limitations in producing a 
needed forecast product. Mellor and Yamada (1974) developed a hierarchy of turbulence 
closure models, ranging from a level 1 model (the smallest equation set with the largest 
number of assumptions) to a level 4 model (the largest equation set with the smallest 
number of assumptions). They tested these equation sets and determined that the level 2 
model performed adequately, with additional improvement going to level 3 that is worth 
the extra computational effort. However the level 4 model, with its 8 additional 
equations, showed little improvement over level 3 for the increased computational effort. 
Mellor and Yamada (1982) further modified their level 3 model, known as the level 2.5 
model, which eliminates the prognostic equation for temperature. The result reduces the 
computational effort with negligible effect to the prediction ability as compared to a level 
3 model. COAMPS™ uses Mellor-Yamada level 2.5 for PBL turbulent mixing and 
diffusion calculations. In this scheme, entrainment at the upper boundary is implicit. 
The vertical profiles of state variables are modified by the TKE profile through the 
turbulent fluxes of heat, moisture and momentum. 

CO AMPS™ uses a turbulence closure model based on the Mellor and Yamada 
(1982) 2.5-level scheme (Hodur, 1997). Following Mellor and Yamada (1982) and Chen 
et al. (2003), the one dimensional TKE budget equation is generally of the form: 


de d —r —— dU —ttSV —— 1 d 

- +—(we) = -wu —~wv ——+ BgwO (wp)-e 

dt dz dz dz p dz 


( 2 . 1 ) 


which is parameterized in COAMPS™ as 
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where e is the TKE, U, V are the horizontal wind speed, K e , Kh and Km are the eddy 
diffusivity constants for TKE, heat and momentum, p is the coefficient of thermal 
expansion, 0 is the potential temperature, and A is an empirical dissipation length scale. 
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In this formulation the prognostic variable, TKE, is determined from the TKE budget 
equation where the production terms are a function of turbulent fluxes. These fluxes are 
diagnosed from the vertical gradient of the scalar quantity times the eddy diffusivity, 
given by 

K,HM,e = ^II.M.e (2.3) 


where S h ,m are polynomial functions of the flux Richardson number, S e is a constant and 
/ is the master length scale. The TKE equation is solved explicitly using a tri-diagonal 
technique to invert the coefficient matrix. Various length scales (e.g. Blackadar 1962, 
Therry and LaCarrere 1983) are available; however the default master length scale 
follows Blackadar (1962) such that for the neutral or convective boundary layer it 
increases with height from a minimum of 5 meters near the surface to a maximum at the 
boundary layer top, and then returns to the minimum value. 


1 = 



(2.4) 


where k is the von Karaman constant, z is the elevation, and a is a stability function based 
on z/L. Some turbulence can be generated at the upper interface in this scheme by wind 
shear and cloud-generated buoyancy flux between grid levels, and this indirectly 
represents the entrainment effects. This approach implicitly represents the entrainment 
mixing across the upper interface at grid resolution with the calculation of turbulent 
fluxes based on the eddy diffusivity and the vertical gradient of the scalar quantity. This 
approach can partially represent the mixing across the interface; however this implicit 
entrainment at vertical resolutions on the order of 100 meters does not accurately 
represent the sharper gradient in the entrainment zone, which is generally sub-grid scale 
of order 10 meters. It is not yet feasible to run operational prediction models at very high 
vertical resolutions, and LES studies indicate that, even at LES resolutions, sub-grid scale 
(SGS) processes may still cause problems with accurately determining the entrainment 
velocity (Moeng 2000; Stevens 2002). Using an adaptive grid scheme to concentrate 
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levels in a diagnosed entrainment zone (Fiedler 2002) showed some advantage in 
increasing vertical resolution without large computational expense, but resulted in no 
improvement in predicting entrainment velocity at resolutions typical for operational 
weather prediction. 

At least for the near future, improvement in mesoscale prediction models of the 
CTBL will require more efficient representation of the sub-grid scale turbulent transport 
near the cloud top. It is hypothesized that this more representative turbulent mixing will 
then feedback on thermodynamic prognostic equations and the cloud microphysics 
scheme to better reflect the observed CTBL in terms of cloud height, thickness, and 
LWP. 

D. PARAMETERIZATION OF ENTRAINMENT VELOCITY 

As discussed in Lilly (2002b) and Stevens et. al (2003b) the PBL 
parameterizations in General Circulation Models (GCM) are usually at relatively coarse 
vertical resolutions and often explicitly solve for boundary layer height and therefore 
must explicitly represent entrainment. The motivation for much of the recent work in 
parameterizing entrainment has been to improve the ability of GCM’s to accurately 
diagnose entrainment velocity with an explicit approach at coarse vertical grid spacing. 
This is required to represent a realistic CTBL and in particular stratocumulus to cumulus 
transition, which is important to the global radiative budget and therefore the 
applicability of GCM’s to climate studies. As discussed in Lilly (2002a, 2002b), the key 
to appropriately parameterizing the entrainment rate is to represent the relative 
contributions to the entrainment velocity (W e ) of turbulent transport due to surface fluxes 
and the buoyancy flux generated by cloud. Comparison to observations and 
representation of the net effective entrainment is complicated in numerical models by the 
horizontal grid resolution in which the real atmosphere may be only partially cloudy, 
cloud top may undulate within a grid box, and the interface may only be tens of meters 
thick or less while the model vertical resolution is typically an order of magnitude larger. 
Stevens (2002) compares several recent entrainment rules proposed by Lock (2000), 

Turton and Nicholls (1987), Konor and Arakawa (2001), Moeng (2000), and Lilly 
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(2002a) and found that entrainment rules differ by as much as a factor of 2 or 3, even 
when based on identical LES data, and that the coarse vertical resolutions found in GCM 
and operational forecast models cannot currently satisfactorily represent the CTBL 
regime. Following Stevens (2002), Equation 2.5 describes the general form of all of 
these rules and they differ primarily in how they account for buoyancy reversal, shear and 
moist thermodynamic effects. Most of the recent rules can be expressed in the form 

We= A(~~) + D (2.5) 

A b 

where A is an efficiency parameter to describe the energy available to drive turbulence, 
part of which works against the stratification, W is the rate of work based on the forcing 
from the cloud-top buoyancy flux, radiative flux and total water flux, Ab is the stability 
condition, or “jump condition” at the inversion, and D describes forcing from diabatic 
processes such as radiative cooling above the inversion, drizzle, and other processes. The 
various rules treat each of these terms slightly differently and also specify different 
empirical constants depending on the observations and LES studies used in their 
development and the environmental conditions where their assumptions and constraints 
remain valid. 


As discussed in Grenier and Bretherton (2001), the entrainment velocity can be 
diagnosed fairly simply by applying a local closure such as the Tumer-Deardorff scheme 
(Deardorff 1972, 1976) given by 


We=A( £. )=A HL 

Ri LA,b 


( 2 . 6 ) 


where Ri is the bulk Richardson number, Ri=U 2 /LAib, L is the velocity scale, U is the 
horizontal wind shear, and A,b is the jump condition across the inversion. In the common 
notation of Equation 2.5, the Tumer-Deardorff closure in Equation 2.6 becomes a simple 
ratio between the turbulent mixing and the inversion strength with D=0 as shown in 
Equation 2.7. 


in 

We = A, —— (2.7) 

1 V A t b 
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A l is the entrainment efficiency, e is the TKE, and /' is the mixing length. The turbulent 
mixing can be defined as the TKE at the grid level at inversion base or as a vertical 
integral of the mixed layer TKE. This parameterization does not account for variations in 
the in-cloud conditions such as droplet size that would affect the rate of evaporative 
cooling due to mixing between inversion air and cloud parcels, the role of wind shear 
above the cloud-top, or the additional radiative or moist-thermodynamic forcing due to 
the presence of cloud. Another approach by Lilly (2002b), given by Equation 2.8 
attempts to better represent the differences between cloudy and clear convective 
boundary layers by including a consideration of the cloud “wetness” as a function of 
LWC and total water and potential temperature jump conditions. This attempts to 
account for buoyancy reversal in a weighted sense and might be expected to accelerate 
the dissipation of the cloud layer, but still sets the diabatic term (D) equal to zero. In this 
formulation 


A = 


, 5 b 

" A b . 

, . 8b ' 

1 + A„ — u 


w = \ A,o+n, Ao+nJW; 
d = o 


( 2 . 8 ) 


where A n is the entrainment efficiency, a semi-empirical constant, 6b is a measure of the 
stability strength in the undulation layer as a function of cloud wetness and resulting 
change in buoyancy due to mixing of saturated and unsaturated parcels, p relates 
entrainment as a proportion of the net forcing and the W term is equated to local 
weighting (r|) of the top-down and bottom-up radiative (F), moisture (q), and moist static 
energy (s) fluxes. Another approach by Lock and MacVean (1999) attempts to include 
diabatic (radiative) effects and buoyancy reversal with 
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where the formulation for W assumes that layer energetics are consistent with what one 
would expect in a uniformly unsaturated layer with an additional forcing term due to 
buoyancy reversal where X* is the amount of free-tropospheric air as related to (1-X) 
parts boundary layer air. The D term addresses the efficiency with which radiative 
forcing and water phase changes drive buoyancy fluctuations. These changes in the W 
and D terms enhance the efficiency with which radiation and entrainment mixing are 
allowed to generate turbulence. 

Based on results from DYCOMS II Flight 1 as reported in Stevens et al. (2003a, 
2003b) and considerations discussed in Stevens (2002) and Grenier and Bretherton 
(2001), our initial parameterization implements the Tumer-Deardorff closure as it can be 
easily expressed within the CO AMPS™ framework. We have also modified the length 
scales and determination of the jump condition as discussed in Appendix B of Grenier 
and Bretherton (2001) to partially account for the evaporative cooling effects caused by 
the mixing fraction between cloudy and clear air. 

The three primary physical processes responsible for the CTBL structure are 
radiation, turbulence, and cloud microphysics. The CTBL has an additional source of 
turbulence not present in a clear convective boundary layer in that radiative flux 
divergence caused by the presence of the cloud results in additional buoyancy flux that 
drives turbulence at the cloud top. If evaporative cooling from mixing entrainment air 
into cloudy parcels is greater than entrainment warming, turbulence forced entrainment 
and can lead to cloud dissipation through buoyancy reversal. This is generally referred to 
as Cloud Top Entrainment Instability (CTEI) (Deardorff 1980, Randall 1980). As 
discussed in Stevens (2002), the relative role various mechanisms such as SST, CTEI, 
and drizzle play in CTBL formation, dissipation, and stratocumulus to cumulus transition 
is still unresolved. One area of active research with a strong effect on the CTBL structure 
is in the proper specification of surface fluxes and turbulent transport in the TKE closure 
(e.g. Miller et al. 2001, Rados et al. 2002). The specification of entrainment fluxes at the 
cloud-top sub-grid scale interface might be expected to have less effect, however this 
mechanism still results in the model being unable to capture boundary layer evolution 
completely as discussed in Grenier and Bretherton (2001) and Lock (2001). 
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III. EVALUATION OF COAMPS™ FOR STRATOCUMULUS 
TOPPED BOUNDARY LAYERS USING OBSERVATIONS 

A. THE DECS FIELD EXPERIMENT AND DECS DATA ANALYSIS 

In order to better understand the role of entrainment in the evolution and structure 
of coastal stratocumulus and the interactions with the complex coastal flow field, our 
approach will consist of three parts: analysis of the model diagnosed boundary layer 
height in comparison with the observed inversion base height, cloud thickness, and 
inversion strength as compared to observations, implementation of a parameterization to 
explicitly represent entrainment and resulting forcing in CO AMPS™, and comparison of 
simulations using this explicit method with simulations using the current formulation and 
with the observed boundary layer structure. This chapter will focus on the evolution of 
the predicted boundary layer field as compared to observations. 

1. Observations - The Dynamics and Evolution of Coastal Stratus (DECS) 
Field Study 

Aircraft, cloud radar, radiometer and rawinsonde observations were collected in 
Central California near Monterey Bay during the Dynamics and Evolution of Coastal 
Stratus (DECS) field study from June 13 to July 22, 1999. Twenty-one daytime research 
flights were flown which made 10 and 100 Hz measurements of marine boundary layer 
turbulence, thermodynamic, cloud microphysical, and aerosol properties using a two- 
engine UV-18A Twin Otter operated by the Center for Interdisciplinary Remote Piloted 
Aircraft Study (CIRPAS) of the Naval Postgraduate School (NPS). The instrumentation 
package, calibration and processing of the turbulence data are described in Kalogiros and 
Wang (2002). Figure 3.1 depicts the research area, surface data sites and two 
representative flight tracks, Figure 3.2 depicts a typical 3-dimensional flight track, and 
Table 3.1 describes the flights analyzed in this study. The first flight had instrument 
problems so the flight numbers in Table 3.1 refer to flights 2 through 21. Researchers 
from the University of Miami, Rosenstiel School of Marine and Atmospheric Sciences 
(RSMAS) operated a fixed, upward looking 94 GHz Doppler cloud radar during part of 
this period that allowed the determination of cloud top height, cloud reflectivity, cloud 
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vertical velocity and drizzle fall velocity. Kollias and Albrecht (2000) provide a 
description of this radar and its processing methods and Fritz et al. (2001) describes the 
use of the radar data to explore drizzle characteristics using a case study from DECS. A 
microwave radiometer operated by Pennsylvania State University was collocated with the 
cloud radar and was used to obtain the integrated water vapor and liquid water path. 
Several rawinsondes per day were launched by NPGS from the radar site on flight days. 
The rawinsondes were modified to measure the lower atmosphere to approximately 4000 
m over about 30 minutes before descending back to the surface and providing an 
additional boundary layer profile during descent. The location of these instruments was 
at approximately 36.699N 121.808W at the City of Marina water treatment site, marked 
“DECS” on Figure 3.1. 

Additional boundary layer observations were archived from the NPS site at Fort 
Ord, approximately 5 km to the northeast. This site is operated by the Naval 
Postgraduate School and includes 915 MHz and 404 MHz wind profilers, a Radio 
Acoustical Sounding System (RASS) to obtain a virtual temperature sounding, a Vaisala 
laser ceilometer, which provides cloud base information since the cloud radar cannot 
make observations below 200 m, and a surface meteorology suite that measures surface 
pressure, wind speed and direction, temperature and dewpoint, and shortwave and 
longwave irradiance. Basic data processing for the 915 MHz profiler is discussed in 
Ralph et al. (1998). 

Surface mesonet data from the Real-time Environmental Information Network 
and Analysis System (REINAS; Nuss et al. 1996) as well as surface METAR, National 
Data Buoy Center (NDBC) and research moored buoy data, and 4km visible and infrared 
imagery and 1km visible geostationary satellite imagery from GOES-10 were also 
archived for the period. 
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Figure 3.1. COAMPS 5km inner grid, terrain, surface data, and representative 
flight tracks for the DECS field campaign. Diamonds indicate surface observation 
sites, buoy data is identified by numbers and the designators Ml-M3. Standard 
synoptic and mesonet surface sites are indicated by three letter call signs. Thick 
solid lines mark two representative flight tracks to the west and south. Flights 
typically consisted of stacked level legs, slant and spiral soundings and sawtooth 
profiles through the lower troposphere. 
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Figure 3.2. A typical three dimensional flight track during DECS. Flights 
included multiple sounding legs through the full extent of the boundary layer, 
“porpoise” or sawtooth soundings providing multiple measurements across the 
inversion, and level leg “stacks” of constant heading and altitude for multiple 
levels below, within, and above the cloud layer. 
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Fli 

Date 

IVinLat. 

IVSnLon 

MaxLat. 

Max Lon. 

Start (UTC) 

Bid (UTC) 

Mu Alt. (rr) 

Max. Alt. (nr) 

1 

618 

36.61 

-122M 

36.71 

-121.66 

1500 

1620 

27 

1500 











2 

619 

36.48 

-12258 

36.99 

-121.75 

1530 

1915 

32 

954 











3 

621 

36.60 

-12257 

36.79 

-121.76 

1430 

1525 

26 

990 











4 

624 

36.66 

-12292 

36.87 

-121.76 

1500 

1900 

21 

967 











5 

625 

36.18 

-12293 

36.74 

-121.60 

1500 

1745 

28 

2403 











6 

628 

35.11 

-12221 

36.80 

-121.04 

1430 

1845 

32 

2140 











7 

629 

35.47 

-12276 

36.78 

-121.15 

1430 

1900 

15 

3266 











8 

702 

34.72 

-12254 

36.71 

-121.00 

1700 

1945 

49 

1414 











9 

703 

34.41 

-121.93 

36.70 

-121.82 

1345 

1815 

25 

1450 











10 

706 

36.63 

-124.39 

36.75 

-121.73 

1500 

1900 

21 

1471 











11 

707 

36.66 

-12276 

36.81 

-121.76 

1540 

1815 

26 

970 











12 

709 

36.66 

-124.38 

36.73 

-121.76 

1500 

1915 

21 

1063 











13 

710 

36.65 

-124.37 

36.73 

-121.76 

1500 

1900 

25 

1410 











14 

714 

36.66 

-124.71 

36.74 

-121.76 

1510 

2000 

15 

1160 











15 

715 

36.65 

-124.70 

36.75 

-121.76 

1500 

1700 

25 

1289 











16 

716 

35.70 

-124.71 

36.74 

-121.76 

1800 

2230 

26 

1357 











17 

717 

36.33 

-123.M 

36.87 

-121.62 

1530 

2015 

28 

1846 











18 

719 

36.67 

-124.71 

36.81 

-121.76 

1450 

2020 

22 

1829 











19 

720 

36.66 

-123.06 

36.74 

-121.76 

1500 

2015 

17 

1827 











20 

722 

36.58 

-125.05 

36.92 

-121.76 

1800 

2200 

29 

1442 


Table 3.1. General information on the DECS flights analyzed in this study. 
Flights included multiple sounding legs through the full extent of the boundary 
layer, “porpoise” soundings providing multiple measurements across the inversion, 
and level leg “stacks” of constant heading and altitude for multiple levels below, 
within, and above the cloud layer. 
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B. COAMPS™ SIMULATIONS FOR DECS CASES 

1. Model Configuration 

In this study, our primary data analysis and visualization software tool was 
MATLAB version 5.3, a program commonly used in the scientific community data 
analysis. Meteorological analysis and forecast simulations were conducted with 
CO AMPS™ version 3.1 as described in Chen et al. (2003). CO AMPS™ was run in a 
mode similar to the current operational forecast model (30 levels), triply nested at 
horizontal resolutions of 45, 15 and 5 km as depicted in Figure 3.3. The vertical grid was 
modified slightly to add resolution in the boundary layer (Table 3.2) and to ensure 
consistent vertical grid spacing across the inversion for ease in data interpretation. 
Several other modifications were also completed and tested in the memory arrays and 
output routines to aid in model result interpretation. These changes create additional 
output variables such as the TKE budget components and radiative and turbulent fluxes. 
Identical comparison runs were done to ensure the code changes designed to produce 
additional output variables did not inadvertently change the memory arrays or model 
results. 

The model was run using the COAMPS™ Ocean Data Assimilation (CODA) in 
order to analyze the ocean sea surface temperature (SST) and ice fields from observations 
at higher resolution suitable for mesoscale forecasts rather than using the analysis fields 
from the global model (NOGAPS). This allows the mesoscale variability in the coastal 
zone from satellite retrievals (approximately 10 km resolution) and ship and buoy data to 
be better represented. Figure 3.4 compares observed SSTs from AVHRR satellite data to 
the model analysis fields. The general features are well depicted including a pronounced 
cooling extending westward and southward of Point Reyes and Point Sur and general 
warming to the south and west overall with coolest temperatures near the coast. The 
model was run in simulation mode for the entire period of the DECS field study, from 
June 13, 1999 to July 22, 1999. This means that the forecast was run to twelve hours and 
then a new data assimilation cycle was run to combine the forecast with observed 
conditions. This is different from the typical operational 48-72 h forecast length in order 
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to better represent the real atmosphere, thus the model tendencies are periodically 
bounded by synoptic observations when comparing the model forecast to observations. 
Observation flights were typically done between 1500 and 2000 UTC corresponding to a 
3 to 8 h model forecast. This period is shorter than the 8-12 hour pre-integration time 
suggested by Cui et al. (1998) for dynamic initialization to ensure all gravity wave 
oscillations caused by the initialization have propagated out of the model domain. In our 
case, however, we decided that capturing the synoptic conditions from the conventional 
observations was more important than ensuring dynamic balance since we are comparing 
model forecasts to specific observations and events. None of the special data from the 
field study was assimilated in the initial analysis for the simulations. Therefore the 
comparisons between the model and the DECS observations are truly independent and 
the mesoscale and marine boundary layer structure in the simulation relies mostly on the 
model physics and large scale forcing since these structures are under-resolved in the 
observations. The COAMPS™ simulations did however include all standard data sets 
used in the operational model such as surface observations, coastal buoys, rawinsondes, 
and satellite data. 


COAMPS nests 



Figure 3.3. 1800 UTC 10m Potential temperature (K) and model nest boundaries 
for the DECS simulations. 
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925 

1100 

i 400 

650 

| 100 

I 940 

750 

300 

! 550 

100 

950 

500 

200 

450 

100 

960 

330 

140 

350 

100 

970 

215 

90 

260 

80 

980 

140 

60 

180 

80 

990 

90 

40 

115 

50 

1000 

55 

30 

65 

50 

1005 

30 

20 

30 

20 

1010 

10 

20 

10 

20 

1015 


Table 3.2. Operational and modified model vertical levels. 
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Figure 3.4. Comparison of CODA Sea Surface Temperature Analysis with 
AVHRR satellite observations for 2 days during the DECS field study. Note the 
general agreement between the model analysis and observations for features such 
as the cooler temperatures offshore of Point Reyes and Point Sur and a warming 
trend to the southwest. 
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2. Modifications to the Calculation of In-cloud Buoyancy Flux 

In initial CO AMPS™ simulations it was noted that forecast LWC and LWP 
values were an order of magnitude too small for an otherwise reasonable thermodynamic 
profile as compared to observations. Investigation into this issue indicated that under¬ 
prediction of low cloud may be attributed to the method used to diagnose of the presence 
of cloud and the applicable thermodynamic variable in calculation of the in-cloud 
buoyancy flux in the boundary layer turbulence scheme. COAMPS™ uses the 
condensation rate to diagnose if a cloud layer is present in calculating the buoyancy 
contribution to turbulence, but this method has the effect of diagnosing a cloud-free 
condition if cloud is dissipating but still present (Shouping Wang 2003, personal 
communication). An incorrect TKE profile therefore results as the fluxes are calculated 
for clear air and the cloud layer thus dissipates too rapidly. Modifications to the 
boundary layer scheme to use a saturation condition vice condensation rate to define the 
cloud layer are planned for incorporation in a future COAMPS™ release. The turbulence 
scheme was modified to identify the cloud layer from the saturation condition and then 
use the liquid potential temperature gradient in-cloud to calculate the stability function 
for the mixing coefficient vice using the virtual potential temperature gradient. This 
produced much better results for the summertime stratus cases in this study as can be 
seen in Figure 3.5. All simulations and COAMPS™ versions discussed in this study 
except for the comparison test in Figure 3.5 use the boundary layer subroutine modified 
as discussed to produce more realistic cloud fields. 
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Figure 3.5. 12 hour forecast of integrated cloud liquid water (kg m' 2 ) for a) 

original CO AMPS™ 3.1 and b) the modified turbulence parameterization. The 
parameterization was modified to 1) identify the cloud layer from the saturation 
condition vice the condensation rate and 2) calculate the stability in-cloud for the 
turbulent mixing coefficient from the liquid potential temperature gradient. 


37 









C. COAMPS™ - DECS COMPARISON 

To better understand the ability of the simulations to accurately forecast the 
boundary layer and cloud structure, the simulation results were compared to aircraft 
observations. Vertical profiles of model output on a one hour interval were extracted 
from the 5 km inner nest for the closest point corresponding to the mean time and 
position for each of 445 aircraft slant sounding legs. The aircraft flights were generally 
between 1500 and 2000 UTC, which corresponds to a 3-8 hour model forecast. The total 
water and buoyancy fluxes from the model were also extracted for each level flight leg. 
Turbulent variances and fluxes can also be calculated using the eddy correlation method 
for the slant path soundings, however the results may be less reliable for sounding legs, 
especially if the aircraft is not maintaining a steady heading and rate of ascent. Mean 
turbulence values from level legs at various altitudes were therefore also calculated and 
compared to slant path and model profiles. Hourly east-west cross-sections for both the 
daytime (1200 UTC) and nighttime (0000 UTC) simulations along the mean flight track 
latitude were also analyzed from the model to assess diurnal development and variability. 
The purpose here is to describe the general CTBL structure and variability in both the 
simulations and observations. This analysis will allow us to determine the performance 
of the current operational forecast model for the period of interest and will be the basis 
for adjustments made to the turbulence parameterization within COAMPS™ via the 
explicit entrainment velocity calculation and fluxes to examine the effect of implicit 
entrainment in the current model formulation. 

The thermodynamic structure of the model atmosphere was found to compare 
well to both the aircraft sounding legs and coastal rawinsondes (not shown). Surface 
temperatures were generally within 1 degree (K) and values in the free atmosphere 
matched well. As shown in Figure 3.6, the east-west cross-sections depict a characteristic 
lowering of the inversion height towards the coast with a much steeper slope within 100 
km of the coast and increasing heights across the coast and the coastal plain until the 
boundary layer intersects the terrain. The observed structure in Figure 3.6.a) does not 
slope as sharply in the inner coastal, which is partly an artifact of the smoothing and 
plotting of limited point data, but still shows the characteristic slope offshore. The wave 
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pattern in the model data and the lowering boundary layer height farther offshore near 
125 W that are not seen in the observations are likely due to a combination of smoothing 
of the observational data and gravity wave propagation from the Sierra Nevada as 
discussed in Cui et al. (1998). The rapidly changing inversion height in the inner coastal 
has a strong diurnal and synoptic characteristic and is associated with terrain-induced 
increased mesoscale subsidence and the MCLLJ. When compared to observations, the 
structure as depicted in Figure 3.6 seems to be well supported, although the slope is too 
steep and the inversion height is too low. Figure 3.7 compares the model inversion 
height based on the vertical 0i gradient to the observations, which shows that this low bias 
is consistent throughout the field study and that there is a mean bias of 200 meters. The 
consistent bias at different forecast times and under different synoptic conditions supports 
our theory that differences between the COAMPS™ simulation and the observed 
structure are not an artifact of a poor initialization or lack of dynamic balance in the 
simulation but rather are caused by inappropriate representation of key physical 
processes. Figure 3.8 is the difference between modeled and observed boundary layer 
height as a function of cross-coast longitude. Although there is more scatter in the inn er 
coastal zone, a dependence on cross-coast distance is not readily apparent. 
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Figure 3.6. Comparison of virtual potential temperature (K) between a) the 
observations and b) the corresponding model forecast. A solid line superimposed 
on the contour plot depicts the flight track for the observations. The COAMPS™ 
simulation is a five hour forecast valid at 1700 UTC, close to the midpoint of the 
aircraft flight. 
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Aircraft Observed PBL Height (m) 


Figure 3.7. Comparison of observed and modeled PBL heights as defined from 
the height of the greatest potential temperature gradient (K m' 1 ) or “jump”. Values 
in the upper left corner are the mean values for the aircraft measurements and 
corresponding COAMPS™ profiles. 
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Figure 3.8. Differences between the observed and COAMPS™ boundary layer 
height normalized by the observed boundary layer height as a function of 
longitude. Boundary layer height is determined from the height of the greatest 
liquid potential temperature gradient (K m' 1 ) or “jump”. (Groupings at 0.2 and 0.5 
correlate to particular flight days and represent the fit between the forecast and the 
synoptic conditions rather than a particular boundary layer regime.) Note the larger 
variability in the inn er coastal zone between 121.5 W and 123.0 W. 
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Figure 3.9 shows a comparison of modeled LWP to the corresponding observed 
cloud field. Qualitatively CO AMPS™ represents variations in the coastal stratus field 
well compared to the satellite observations. Climatologically, the marine stratus is 
relatively persistent offshore with thinning/ dissipation and a more cellular structure 
developing frequently during the day. A thin finger of fog and stratus is often observed, 
which persists in a narrow zone within 30 km of the coast, especially downwind from 
prominent points such as Pt. Arena and Pt. Sur. A characteristic oval shaped cloud 
thinning associated with the onset of the LLCJ is often observed downstream from Cape 
Mendocino and occasionally downstream of less prominent points such as Point Reyes, 
Pescadero Point, and Point Sur. This lowering of the boundary layer and cloud thinning 
often persists for 2-3 days and the diurnal signal is superimposed on this temporally 
synoptic and spatially mesoscale variability. At the coast over land, the cloud patches 
that persist during the day often deepen and develop in late afternoon and become 
indistinguishable from developing nocturnal valley fog, which generally dissipates by 
mid-morning. The model depicts these major observed cloud features quite well 
including the persistent marine stratus, diurnal variation in cloud fraction and thickness 
west of about 123 N, inner coastal clearing associated with the coastal jet and sea-breeze 
forcing and more persistent cloud upwind of headlands, and nocturnal valley fog clearing 
by midday. Figure 3.10, CO AMPS™ is an East-West cross-section of cloud liquid water 
mixing ratio at 37.1 N along the aircraft flight path. A lowering of the cloud top and base 
in the narrow inner coastal zone is seen in both the observations and model output which 
matches our expectation from considerations of lower SSTs and increased mesoscale 
subsidence due to flow interaction with the terrain. The model frequently produces fog in 
this zone (not shown), however observations indicate that, although the cloud top and 
base in this region is generally lower than elsewhere, it is usually not saturated all the 
way to the surface, especially during daytime. 

Figure 3.11 compares the observed and modeled cloud top height. The 200 m low 
height bias seen in the liquid potential temperature gradient is also seen here which 
further supports our conclusion that the boundary layer height is too low in COAMPS™ 
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in the CTBL regime. Figure 3.12 is a comparison between the observed maximum cloud 
liquid water mixing ratio, which occurs at cloud top, and the modeled values. 
COAMPS™ may be expected to have slightly lower peak values than the observations as 
a result of representation of a sharp gradient on a coarse discrete grid. Instead, this 
comparison shows that the maximum liquid water mixing ratio in COAMPS™ is 
consistently about 25% too high. The low cloud top bias and high cloud LWC bias have 
important implications to the radiation budget, optical depth, and surface heating by solar 
warming in the model and are a major focus of this study. Figure 3.13 is a comparison of 
the cloud base to ceilometer data for a representative two day period. The ceilometer is 
sufficiently sensitive to scattering by liquid water and aerosols to detect liquid water 
mixing ratios to about 0.02 g kg" 1 . Based on the ceilometer and aircraft data, both the 
cloud top and base in CO AMPS™ are about 200 meters too low. The cloud also lifts and 
dissipates 1-3 hours later than indicated from observations and reforms 1-3 hours sooner, 
which is consistent with over-prediction of LWC. 
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Figure 3.9. Vertically integrated liquid water mixing ratio (g kg" 1 ) (LWP) from 
COAMPS™ and corresponding visible 1 KM satellite imagery. The model 
qualitatively matches the observations for major coastal cloud features including: i) 
persistent marine stratus with a diurnal change in cloud fraction, ii) coastal clearing 
with persistent patches upwind from points, and iii) coastal stratus and inland 
valley fog with a strong diurnal character. 
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Figure 3.10. Comparison of cloud liquid water mixing ratio (g kg' 1 ) between the 
model and the observations along and East-West cross-section at 37.1 N, which 
corresponds with the aircraft flight track. A solid line superimposed on the contour 
plot depicts the flight track for the observations. The CO AMPS™ simulation is a 
five-hour forecast valid at 1700 UTC, close to the midpoint of the aircraft flight. 
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Figure 3.11. Scatter plot of observed and modeled maximum liquid water mixing 
ratio. 
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Figure 3.12. Comparison of cloud top height (m) between the model and the 
observations. In the case of the discrete model grid, cloud top height is taken as 
the level of maximum cloud liquid water mixing ratio. In the observations, the 
vertical distance between this value and a zero value is tens of meters while the 
model may represent the cloud top across 2-3 grid levels at 100 m vertical 
resolution. 
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COAMPS™ Simulation and NPS Ceilometer Observations 06/23/1999-06/25/1999 



Figure 3.13. Time series plot of observed ceilings (line data) at the NPS Ft. Ord 
ceilometer site and modeled liquid water mixing ratio (g kg' 1 ) (fdled contours). 
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We also examined the predicted wind field relative to the MLLCJ and compared 
it with available observations. The horizontal structure in the simulated MLLCJ depicted 
in Figure 3.14 matches the conceptual model of terrain induced variability well, although 
the horizontal velocity patterns cannot be verified from observations. The presence of the 
jet is often difficult to detect in surface observations in that it reaches its maximum 
magnitude at the boundary layer top, therefore the reflection of the LLCJ in buoy 
observations depends on momentum mixing in the boundary layer. Detection in aircraft 
observations requires that the available data are sufficiently dense and correctly located in 
time and space to capture the complete structure. The horizontal variability in the LWP 
depicted in Figure 3.14 correlates well with the conceptual model of the cloud patterns 
and jet structures arising from the boundary layer response to the coastal flow and 
typically observed cloud patterns seen in satellite data under this regime (not shown). The 
vertical structure of increasing winds at the boundary layer top along an east west sloping 
inversion matches the conceptual model and the aircraft observations well as can be seen 
in Figure 3.15. The aircraft slant path soundings provide a unique opportunity to directly 
study the vertical structure of the jet at varying distances offshore. 

In comparing the observed and modeled vertical wind structure in jet regions it 
was noted that the classic LLCJ velocity profile of a distinct velocity maximum at the 
inversion and decreasing wind speeds towards the surface was rarely observed in cloudy 
cases in the observations. Instead the aircraft soundings often depicted uniform stronger 
winds within and below the cloud and decreasing winds above the cloud top. 
Occasionally wind speed would increase again to a secondary maximum well above the 
boundary layer top. In the no-cloud cases the observations generally matched the 
expected vertical wind speed structure from LLCJ theory of increasing winds in the free 
atmosphere towards the surface to a maximum at the boundary layer top and decreasing 
winds in the boundary layer. On the other hand, the model generally produced a relative 
velocity maximum at the boundary layer top in both cloudy and clear cases as depicted in 
Figure 3.16. 

Figure 3.17 is an example from one flight of the observed and modeled fluxes for 
a series of “stacked” level legs. Because the turbulent fluxes derived from 10 Hz slant 
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profiles of aircraft data using the eddy correlation method have the potential for large 
uncertainty, especially if the aircraft maneuvers or does not maintain a sufficiently slow 
and constant ascent rate, discussion of the observed fluxes is based a subset of vertical 
profiles with a long ascent path and constant heading and averaged flux values calculated 
from stacks of long, horizontal legs. The observations generally showed small surface 
buoyancy fluxes on the order of 10-15 W m at the surface and decreasing to a minimum, 
or even slightly negative buoyancy at the top of the sub-cloud mixed layer. Buoyancy 
flux then increased in the cloud layer to a maximum near cloud top approximately equal 
or larger than the surface buoyancy flux. There was a great deal of scatter in the cloud 
top flux measurements with large positive and negative fluxes occurring. This general 
structure is well correlated between the individual soundings and level leg averages. 
Figure 3.18 and 3.19 depict a typical sounding leg and corresponding model sounding for 
a cloudy case and a no-cloud case respectively. Note the much closer agreement between 
COAMPS™ and observations for the no-cloud in Figure 3.19. Frequently the simulation 
and the observations were in reasonably good agreement in the sub-cloud layer, however 
the model cloud layer was too low, making it difficult to distinguish between surface and 
cloud forcing. In the entrainment zone at and just above the inversion, a narrow region of 
strong negative buoyancy flux was often observed, especially when vertical wind shear 
was large. This was also represented in the simulations, although less frequently and 
with a smaller magnitude. Mixing between the dry inversion air and the cloud can lead to 
evaporative cooling and positive buoyancy flux through buoyancy reversal. Since, the 
turbulence parameterization implicitly represents the entrainment, the tendency of the 
model to under-predict the magnitude of the large positive fluxes in-cloud and negative 
fluxes in the entrainment zone supports our hypothesis that the implicit approach cannot 
adequately represent the entrainment velocity across a narrow sub-grid scale entrainment 
zone. 

In summary, the model captured the structure and evolution of the boundary layer 
reasonably well, with the exception of a consistent high bias in cloud liquid water and 
low bias in inversion strength and height. These biases in the inversion strength and 
height were seen more strongly in cloud-topped cases than in clear boundary layers and 
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were seen more strongly in the potential temperature profile than in the water vapor 
profile. On average, the height bias was 200 meters, approximately 25% of the total 
boundary layer height, but varied from less than 100 meters to more than 600 meters. 
These trends were confirmed with over 400 individual sounding legs, 80 level legs and 
20 flights in varying synoptic conditions. This strengthens our hypothesis that sub-grid 
cloud top processes are a likely the source of this bias, which is investigated in Chapter 
V. 

This study was designed to examine the performance of COAMPS™ at the 
relatively coarse operational vertical spacing of 100m. However, because our reasoning 
for investigating an explicit entrainment approach partly relies on the under-resolution of 
the entrainment zone in the implicit approach, sensitivity tests were conducted to 
examine the model dependency on vertical grid resolution. COAMPS™ simulations were 
run with 40 vertical levels (vice 30) with 50 m resolution to a height of 1000 m in SCM 
and full three dimensional versions. The boundary layer height grew more rapidly in the 
40 level SCM and the cloud liquid water mixing ratio was 10% smaller (not shown). The 
three dimensional tests showed similar results with a sharper cloud liquid water gradient 
as seen in Figure 3.20 that qualitatively matched observed typical CTBL structure better, 
and higher boundary layer heights by 50-100m as seen in Figure 3.21. As will be 
discussed in Chapter VI, the magnitude and direction of these changes are similar to the 
effects of explicitly parameterizing entrainment. The computational expense is much 
larger however, with the EEP having negligible effect on run time while the increased 
vertical resolution increases the wall clock time of a simulation by 50-75%, which in our 
configuration results in a 12 hour forecast taking 16 hours to complete vice 9 hours for 
the control or EEP simulations. Developing a framework in COAMPS™ for EEP is also 
useful for studying the underlying physical mechanisms involved in cloud-top 
entrainment aside from the motivation of improving the operational forecast model. 
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COAMPS™ Simulation Time 1999070617 



Figure 3.14. Predicted Liquid Water Path (LWP) (shaded) (kg m' 2 ) and wind 
speed at the inversion (dashed bold lines) (m s' 1 ). The COAMPS™ simulation is a 
five-hour forecast valid at 1700 UTC, close to the midpoint of the aircraft flight. 
Note that in the COAMPS™ simulation the relative minima in LWP correspond 
with maxima in the wind speed. This horizontal structure is difficult to depict in 
observations due to our inability to sufficiently resolve the large scale three 
dimensional structure in observational data, however the model depiction matches 
well with satellite observations of cloud patterns and point observations of wind 
velocity. 
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Figure 3.15: Comparison of wind speed (m s" 1 ) between the model and the 
observations. A solid line superimposed on the contour plot depicts the flight track 
for the observations. The COAMPS™ simulation is a five-hour forecast valid at 
1700 UTC, close to the midpoint of the aircraft flight. 
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Forecast Time: 1999071017 UTC Observation Time: 17.02 Sounding 23 








Wind Direction (deg) Location and PBL Height (m) 


Figure 3.16. In cloudy cases the vertical profile of horizontal velocity (lower left 
panel) was often uniform with height in the boundary layer and dropped off rapidly 
above cloud top. COAMPS™ simulations frequently depicted a velocity 
maximum just above the boundary layer top, which matches the clear profile more 
closely. 
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Observed vs. Model Level Legs 19990717 16-17.7UTC36.7-122.59 



Net Heating Rate (K day 1 ) 







TKE (m 2 s -2 ) 


Figure 3.17. Typical profiles of observed values for aircraft level legs. Small 
scattered dots are individual data points and circles are mean values. 
Corresponding model profiles are solid lines marked with a diamond at each grid 
level. Note the in-cloud buoyancy flux maximum and sub-cloud minimum are 
under-resolved by CO AMPS™. Net radiative heating rate (upper left panel) is not 
available from the observations. 
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Forecast Time: 1999071618 UTC Observation Time: 18.13 Sounding 2 









Location and PBL Height (m) 


Figure 3.18. Typical model vertical profile (diamonds) and sounding profile (dark 
thick line) for a cloudy case with a minimal coastal jet. 
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Forecast Time: 1999070717 UTC Observation Time: 17.28 Sounding 13 








Figure 3.19. Typical model vertical profile (diamonds) and sounding profile (dark 
thick line) for a no-cloud case with a strong coastal jet. 
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Figure 3.20. East-West cross-section at 37.1 N for: a) 30-level CO AMPS™ at 
approximately 100 m vertical spacing in the boundary layer and b) 40-level 
CO AMPS™ at approximately 50 m vertical spacing in the boundary layer. 
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Figure 3.21. Comparison of CTBL structure for 40-level and 30-level simulations 
as discussed in the text. 


60 













































































































IV. DEFINING THE BOUNDARY LAYER HEIGHT 


A. OBJECTIVE 

One of the most important aspects in being able to implement an explicit 
entrainment parameterization in COAMPS™ is to accurately identify the boundary layer 
top at which entrainment fluxes will be applied. This is challenging for a mesoscale 
model since the theoretical criteria for dynamic stability that use either turbulence fluxes 
(flux Richardson number) or vertical gradients (gradient Richardson number) at the 
interface generally are not applicable for the coarse vertical resolution used in mesoscale 
models. We will first examine the diagnostic boundary layer height using the existing 
method in COAMPS™ against those observed in DECS. Efforts are then made to search 
for alternative approaches to improve the diagnostic boundary layer height. This chapter 
focuses on better identification of the boundary layer top as a diagnostic quantity; 
however improvement of the prediction of the boundary layer properties and thus the 
boundary layer height can only come from modifications to the model numerics or 
physics. Modifying the model numerics with regard to entrainment parameterization is 
discussed in the chapters to follow. 

As described in Chen et al. (2003) the diagnostic boundary layer depth is based on 
the lowest level at which the bulk Richardson number exceeds a critical value, which is 
set at 0.5 in the model. The theoretical basis for this is that when the flux Richardson 
Number is greater than a critical value, generally taken to be 0.25, the flow becomes non- 
turbulent. When the equations are cast in a finite difference form, the appropriate critical 
value at discrete grid levels becomes less certain. Previous research on MABL structure 
and evolution during the Coupled Boundary Layer Air-Sea Transfer (CBLAST) field 
study in Cape Cod, MA (Qing Wang, 2003, personal communication) indicated large 
discrepancies between model diagnosed boundary layer height based on the Richardson 
number (turbulence profile) and the inversion base height indicated by manual estimation 
from model sounding profiles. By examining a variety of boundary layer and inversion 
conditions and comparing modeled structure to aircraft sounding data, we refine the 
diagnostic method to ensure accurate representation in the model. 
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B. COAMPS™ DIAGNOSTIC BOUNDARY LAYER HEIGHT 

The boundary layer top is the interface between the turbulent boundary layer 
below and non-turbulent tropospheric air above. As a result, the boundary layer top is 
characterized by a significant increase in Richardson number. Theoretically, this 
interface is best identified using the flux Richardson number: 


R f = 


e 


— dU — 

wu -h W V 


dz 


dT 

dz 


(4.1) 


where the numerator and the denominator are the buoyancy and shear production terms 
of the turbulent kinetic energy (TKE) budget equation, g is gravity, 0 V is the virtual 

potential temperature, U and V are the mean horizontal wind components, and u’, v’, w’ 
are the velocity perturbations. The flux Richardson number defines the dynamic stability 
of the flow. Under stable stratification, turbulence will cease if the buoyancy term is 
larger than the shear production term, specifically when R f >1. Thus a critical flux 

Richardson number of 1 can be used to define the dynamic stability of the flow. 
Alternatively, we can define a gradient Richardson number: 


R, = 


g_dd_ 
0 dz 


fduY 

v dz J 


+ 


dr 

dz , 


(4.2) 


where the gradient is used instead of the fluxes. The critical gradient Richardson number, 
R ic , was deduced to be !4 for inviscid stability by Miles (1961). Measurements suggests 

that the range of R ic is R ic ~ 0.2 ~ 0.25 (Garratt, 1992). In a numerical model with finite 

vertical resolution, the gradient Richardson number is replaced by the bulk Richardson 
number defined as: 


62 



= Ad x Az 
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(4.3) 


* (At/) 2 + (AF) 2 

Because of the bulk approximation, the critical bulk Richardson number, Rb c , is often not 
the same as Rj C . In fact, when the grid resolution is coarse (i.e., large A z) sub-critical 
Richardson number values are almost certainly present in thinner layers within the 
vertical grid interval, implying that there are sub-grid scale fluxes supported by the 
smaller-scale turbulence. This effect is allowed for by increasing the magnitude of the 
critical Richardson number as is currently done in COAMPS™ (Chen et al., 2003) where 
the critical bulk Richardson number is set at 0.5. In addition, COAMPS™ interpolates 
linearly between grid levels to estimate the exact height where the Richardson number 
exceeds the critical value. Thus, the diagnosed boundary layer height may not be at the 
grid levels. 

As was discussed in Chapter III, comparison between the CO AMPS™ diagnosed 
boundary layer height and observations indicates significant under-prediction of the 
boundary layer height. We theorize that this discrepancy is likely because either: 1) the 
boundary layer height was diagnosed inappropriately from the model output or 2) the 
boundary layer properties were predicted inappropriately. The work in this chapter 
intends to identify and improve the error caused by an inappropriate boundary layer top 
diagnostic scheme. 

To find a better diagnostic indicator for boundary layer height, one needs a 
variable independent from the bulk Richardson number. For the stratocumulus-topped 
boundary layer, a natural choice is the cloud top. Since liquid water generally increases 
with height in stratiform marine clouds and rapidly drops to zero just above the maximum 
value (Albrecht et al. 1988), the cloud top can be considered as the height where cloud 
liquid water reaches a maximum. Figure 4.1 shows a comparison between the height of 
maximum liquid water content (cloud top height) and the boundary layer height from all 
DECS soundings made by the Twin Otter. Here the boundary layer height was selected 
manually from the height at which the largest gradients in total water mixing ratio and 
liquid water potential temperature are co-located. As can be seen in Figure 4.1, the two 
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variables are nearly the same with a mean difference of about 12 m. Thus the 
observations show that the cloud top height is indeed a good indicator for boundary layer 
height in this regime. Since were are primarily concerned with the CTBL in this study, 
this is a useful indicator for boundary layer height, while another approach such as the 
stability criteria or gradient strength would be more appropriate for the clear case. 


C. COMPARISON OF OBSERVED AND SIMULATED BOUNDARY LAYER 

HEIGHT 

Using cloud top height as a measure, boundary layer height was diagnosed using 
several different indicators including the Richardson number criteria and was compared 
to a variety of observed and modeled thermodynamic variables including the cloud top 
height, and the vertical gradients of liquid potential temperature gradient (d0i/dz) and 
water vapor mixing ratio (dq t /dz) based on the fact that largest gradients of these are at 
the boundary layer top in the typical marine CTBL. We refer to these methods as the 
inversion strength methods and the cloud top height method for consistency. The 
comparison shown in Figure 4.2 was done between the model and the observations for 
450 model forecast soundings valid at the time and location of the sounding flight legs 
for the twenty DECS flights. The Richardson number-based method interpolates between 
grid points to estimate the exact height where the Richardson number exceeds the critical 
value while the inversion strength method is at the height of the closest grid point. 
Consequently, overlapping/clustering of data points at grid levels occurs so some of the 
spread is due to the revised methods being diagnosed on the discrete grid levels and the 
current method being interpolated between levels. However, even with taking the 
variance due to interpolated versus discrete grid level into account, the thermodynamic 
gradient methods show improvement. Vertical grid spacing is about 100 m as discussed 
in Chapter III. 

Testing the sensitivity of the resultant boundary layer height on choices of the 
diagnostic criteria was done for each of the methods discussed above. Results shown in 
Figures 4.2b and 4.2c are considered the best results for the inversion strength methods, 
where an empirical value of 0.02 Km' 1 were chosen as the criterion for 6 , gradient and 
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an empirical value of -0.0065 g kg' 1 m' 1 was chosen for the q t gradient. These critical 
gradient values are empirically derived for this dataset and would likely be sensitive to 
the strength of the subsidence and surface forcing when the CTBL is at approximate 
equilibrium. 

The comparison between the cloud-top height and the other boundary layer 
diagnostics schemes suggests that the 0 t gradient method diagnoses boundary layer 

heights that are closest to the cloud top height with a difference of less than 20 m (Figure 
4.2b). The Richardson number based boundary layer height is in general lower than the 
cloud top height by about 128 m (Figure 4.2a), while the boundary layer height based on 
total water gradient over-estimated the boundary layer height by about 80 m and had the 
largest scatter (Figure 4.2c). In the observations, excellent correlation was found 
between the maximum cloud liquid water mixing ratio and the inversion base, as 
expected from our conceptual model of the general structure of stratocumulus-topped 
boundary layers (e.g. Brost et al. 1982a and 1982b, Albrecht et al. 1988, Paluch and 
Lenschow 1991). 

In comparing the observed inversion height to the modeled height, it was also 
found that, in addition to the low bias of 200 meters discussed previously, under certain 
conditions the Richardson number based method diagnoses a height 600-700 meters 
lower than observed. Further investigation determined that in these cases, a weakly 
stratified sub-cloud layer existed that met the criteria for non-turbulent flow. An example 
of the temperature and moisture structure in these cases is given in Figure 4.3. This 
structure could occur in a decoupling case, when the turbulence in the cloud layer is 
maintained by in-cloud forcing and non-surface generated mixing. While this height 
meets the definition of being the first non-turbulent layer from the surface, and may be 
useful in other applications such as determining electromagnetic propagation conditions, 
our interest is in reliably diagnosing the main capping inversion at cloud top. Using the 
inversion strength method would eliminate the false boundary layer definition in 
decoupling conditions. 
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D. AN ALTERNATIVE METHOD FOR DIAGNOSING THE MODEL 

BOUNDARY LAYER HEIGHT 

Based on these results, an alternative method was formulated and implemented in 
COAMPS™ using the gradient method, which calculates the liquid potential temperature 
gradient (0i) to diagnose the boundary layer height and then bases the inversion height on 
the lowest level at which d0i /dz exceeds a critical value. This method seems to work 
well for the cloudy boundary layers where strong temperature gradient is found at the top 
of the boundary layer. However, it may be problematic in clear conditions or when the 
boundary layer top has a weak temperature gradient. Since the Richardson number-based 
boundary layer height works rather well for the clear boundary layer, for implementation 
in the 3-D COAMPS™, we used a hybrid approach that adopts the gradient for the 

cloud-topped boundary layer and the Richardson number for the clear boundary layer. 
As an alternative, we simply used the cloud top height, defined as the level with 
maximum cloud liquid water, as the boundary layer height for the cloudy condition. This 
final approach was used for implementing explicit entrainment parameterization. 
Although the 0i gradient method worked well as a diagnostic in the control runs, in single 
column testing of the explicit entrainment parameterization we found that any method 
which uses either a turbulence-based or thermodynamic gradient stability criteria would 
oscillate between grid points when the boundary layer was changing in height. This 
oscillation created stability problems and unrealistic entrainment velocities in the forecast 
model so the boundary layer height diagnostic was modified to determine the inversion 
from cloud top height and a relaxed critical Richardson number-based criterion above the 
level of maximum liquid water content. A logical, consistent, and stable method for 
determining boundary layer height on a discrete grid for the time-evolving boundary 
layer still needs further research. Furthermore, the methods used here are particular to our 
regime of interest and focus on the marine CTBL. A method that works equally well in 
identifying the boundary layer height for the CTBL, the dry convective boundary layer, 
and stable boundary layers over land and water is required, although the current 
Richardson-number based method is expected to work better in non-cloudy boundary 
layers without sub-cloud weakly stable layers. 
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DECS 



Figure 4.1. Comparison of the observed cloud top height and the observed 
inversion height based on the vertical gradient of virtual potential temperature. 
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COAMPS COAMPS 




COAMPS cloud top height (m) COAMPS cloud top height (m) 


Figure 4.2. Comparison of a) the current method of determining the boundary 
layer height based on the lowest level at which the bulk Richardson number 
exceeds a critical value and revised methods based on vertical gradients of b) 
liquid potential temperature and c) total water mixing ratio exceeding an 
empirically derived critical value. In the case of weak inversions that do not meet 
the criteria, the new method defaults to the original method. Note that the liquid 
potential temperature gradient shows the best match with cloud top height and does 
not mis-diagnose the sub-cloud inversion as the boundary layer height. 
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July 14, 1999 180 UTC Sounding 14 




Figure 4.3. Example of mis-diagnosis of the boundary layer height using the 
current method in CO AMPS™. Diamond indicates the diagnosed height from the 
revised method based on the 0i or q v gradient. Plus sign is the result from the 
original method. 
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V. EXPLICIT ENTRAINMENT PARAMETERIZATION 


As discussed in Chapter II, the current COAMPS™ boundary layer 
parameterization represents the turbulent entrainment at the boundary layer top through 
the K closure, which we refer to as the implicit method. Such an approach has 
advantages because it is simple to implement since it requires no special treatment of the 
boundary layer top and the inversion thermodynamic structure, and it is a TKE-forced, 
physically based parameterization. However, it may not be able to deal with entrainment 
effectively without changing the mixing within the boundary layer. An effective 
treatment for the entrainment process is to implement an explicit entrainment 
parameterization (EEP). Implementing EEP provides us a framework that allows us to 
focus on the boundary layer top in order to realistically represent the entrainment process. 
In addition, EEP allows us to more easily test new parameterizations for entrainment 
velocity as they emerge. 

There are three components to implementing an EEP: determining the altitude at 
which to apply entrainment flux; specifying AX, X being a predicted mean quantity; and 
incorporating the entrainment fluxes into COAMPS™. In Chapter IV, we discussed 
various methods for correctly identifying the boundary layer height. In this chapter, we 
will focus on defining the boundary layer top jump conditions and implementation of 
EEP into COAMPS™. Results from testing the EEP implementation in single colu mn 
model (SCM) COAMPS™ will be presented to illustrate the role of EEP in modifying 
the boundary layer cloud and thermodynamics. 


A. EXPLICIT PARAMETERIZATION CONCEPT AND FORMULATION 

As discussed in Chapter II, representing entrainment flux at the CTBL top 
requires specification of both the entrainment velocity (W e ) and the inversion jump 
condition (8X). The latter is associated with changes in the quantity (X) across the 
inversion base. Figure 5.1 gives an example of the vertical variation of potential 
temperature and water vapor mixing ratio typically found across the inversion for the 
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CTBL. In this example, the observed change in virtual potential temperature is about 12 
K in a layer of about 70 m in depth. The decrease in total water is about 4 g kg' 1 within 
the same depth. These sharp variations justify the assumption of a discontinuity or 
“jump” condition usually implemented in one-dimensional mixed layer models. In 
reality, the inversion height undulates horizontally due to surface forcing inhomogeneity 
and mesoscale convective circulations (Atkinson and Zhang 1996), which may result in 
discrepancies between the model representation of the jump height as compared to point 
observations. Evidentially, such sharp gradients cannot be represented in operational 
mesoscale forecast models such as COAMPS™ which typically use a vertical spacing of 
100 m or larger. To properly represent a jump across a 50 m entrainment zone, one would 
expect to need vertical resolutions on the order of 10-25m (2-4 AZ). 

Because of the presence of the sharp inversion, the representation of entrainment 
fluxes becomes very sensitive to the vertical resolution of the model. Lenderink and 
Holtslag (2000) demonstrated that simple turbulence closure models such as a 1.5 order 
prognostic-TKE diagnostic-length scale (E-l scheme) formulation can adequately 
represent entrainment at vertical spatial scales of 25 meters and temporal scales of 1 
second, but degenerate rapidly when approaching the relatively coarse resolutions 
currently used in operational mesoscale prediction models of 100m and 120s. 
Consequently, special considerations are needed in representing the inversion jump 
condition in addition to using appropriate entrainment parameterizations. Lenderink and 
Holtslag propose representing the sub-grid jump condition by looking at the total energy, 
which is conserved, but requires taking into account the radiative flux divergence in 
addition to the thermodynamic fluxes. 

Grenier and Bretherton (2001) proposed three specific methods employing 
conservation of total energy for representing entrainment at the inversion height at 
operational grid resolutions. As depicted schematically in Figure 5.2, these are the 
“prognostic inversion” using MABL depth as a prognostic variable, the “reconstructed 
inversion” using grid level thermodynamic variables to reconstruct the MABL height and 
capping inversion strength, and the “restricted inversion” forcing the inversion to lie on a 
flux level of the host model grid. The first method is not appropriate for models with 
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fixed vertical grid levels such as CO AMPS™ and the third method, although the simplest 
to implement, did not perform well for the CTBL case using SCM tests of a 1.5 order 
turbulence closure model (TCM). We initially implemented the reconstructed inversion 
method in SCM testing in COAMPS™, however, we were unable to produce a 
numerically stable three dimensional model, so we reverted to the simpler restricted 
method during initial development. 

Both of these methods have been tested by Grenier and Bretherton (2001) with 
generally good results, although the reconstructed inversion method was found to be 
highly sensitive to the vertical advection scheme. In Grenier and Bretherton (2001), the 
reconstructed inversion method combined a 1.5 order TCM with an explicit entrainment 
parameterization by using the TCM to determine the fluxes at all levels below the grid 
level (M), which is the closest model level at or below the inversion height (i). For the 
next highest flux level (M+l/2), the entrainment velocity (W e ) is first calculated at the 
inversion using the Tumer-Deardorff entrainment parameterization with a local closure 
assumption as discussed in Chapter II and depicted in Equation 2.5, repeated here as 
Equation 5.1 for clarity. 

3/ 2 

We = A, —— (5.1) 

1 V\b 

This entrainment rate is then related to the entrainment flux through Equation 5.2 below 
and extrapolated to the model flux level. The fluxes at level M+3/2 are assumed to be 
zero by definition. The term “reconstructed inversion” refers to the method of 
determining the appropriate values for the state variables at the inversion based on values 
at discrete grid levels. In this case a height weighted average based on the diagnosed 
inversion height (i) of the values at X M and X M +i is used for any scalar (X). In the 
restricted inversion method, the jump condition and the radiative flux and thermodynamic 
values for the entrainment parameterization are simply determined from the values at grid 
level. 

As reported in Lock (2001), Lock et al. (2000), and Martin et al. (2000), the 
reconstructed inversion method has been implemented in the United Kingdom 
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Meteorological Office (UKMO) Unified Model (UKUM) with good results on the GCM 
scale in depiction of StCu to Cu transition and mixed results on the mesoscale in 
depiction of frontal and orographically forced stratus. The ABL parameterization in the 
UKUM is based on non-local TKE closure following Holtslag and Boville (1993) where 
the TKE profile is prescribed and the magnitude is diagnosed based on the surface 
forcing. The original scheme is modified by an additional profile to account for mixing 
driven by cloud top radiative cooling following the work reported in Van Meijgaard and 
van Ulden (1998). This entrainment scheme is modified as discussed in Lock (2001) to 
adjust for unphysical entrainment effects due solely to the numerical scheme resulting 
from calculating entrainment on separate grid levels from the model radiation, cloud 
microphysics and divergence. This effect, termed “numerical entrainment” is discussed in 
more detail in Lenderink and Holtslag (2000). Using a 1.5 order TCM in SCM tests as 
compared to LES results, Lenderink and Holtslag (2000) found that, at the coarse vertical 
resolutions found in operational prediction models, under-prediction of turbulence was 
compensated for by the numerics of representing large-scale subsidence and clouds in 
Eulerian grid boxes so that the resulting entrainment rates looked reasonable. This is 
highly undesirable in an operational setting and should be taken into account in explicit 
entrainment parameterizations. In our formulation we account for cloud top radiative 
cooling and evaporative enhancement but the sensitivity of the resulting boundary layer 
height in the mesoscale prediction model to the vertical advection scheme under 
subsidence discussed in Lenderink and Holtslag (2000) requires further research. 
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Forecast Time: 1999071016 UTC Observation Time: 15.58 Sounding 4 




Figure 5.1. Example of observed (solid line) and modeled (diamond marker) 
profiles for a) water vapor mixing ratio (g kg' 1 ) and b) virtual potential temperature 
(K) for the stratocumulus topped boundary layer. 
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Figure 5.2. Sketch of the three approaches described in the text for treating the 
inversion in a numerical model. Dashed lines indicate flux levels, x’s indicate grid 
levels. In (b), O denotes the profde of another arbitrary conserved scalar X. From 
Grenier and Bretherton (2001). 
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B. TESTING EEP IN SCM COAMPS™ 

Following Grenier and Bretherton (2001), an explicit entrainment 
parameterization was coded and tested in COAMPS™. Once the framework was 
developed, various methods for identifying the boundary layer height, parameterizing the 
entrainment velocity, and introducing the resultant fluxes back into the forecast model 
were tested. Our general approach was to diagnose the jumps in vertical gradients of the 
virtual liquid potential temperature gradient (0 v i) and total water (q t ), which are conserved 
for moist adiabatic processes, and calculate the entrainment velocity from Equation 5.1. 
The entrainment flux for any scalar (X) can then be calculated from Equation 5.2. At the 
boundary layer top, the relationship in Equation 5.3 can be applied and once the 
entrainment fluxes for heat and total water are determined, a new eddy diffusivity 
constant (K H inv) can calculated by inverting the flux-profile relationship as in Equation 
5.4. CO AMPS™ does not explicitly use the turbulent fluxes to predict the mean variables 
but rather calculates the eddy diffusivity (Kh, Km, K e ) based on the predicted TKE from 
the budget equation. The resulting TKE, mixing lengths, vertical gradients of the scalar 
quantities, and surface fluxes are then passed to a linear solver which uses a tri-diagonal 
matrix to modify the vertical profile of the mean scalar quantities. Therefore, the 
entrainment diffusivity constant can then be used in place of the eddy diffusivity constant 
from the main TCM for both the prognostic equation in the TKE budget and in the scalar 
mixing routines for the next timestep. 

Once the entrainment velocity is diagnosed, the entrainment flux for any 
scalar (X) is then determined from: 

w ' X ' = -W e A,X (5.2) 

The effects of the explicit entrainment parameterization are introduced into the 
TCM in COAMPS™ by using the relation: 

_ AV 

w'X' inv =-K h inv — (5.3) 

Az 

where K H inv is the eddy diffusivity at the boundary layer top. Combining Equation 5.2 
and 5.3, we have: 


77 



(5.4) 


K = W Az 

^11 INV '' e^INV 


This is then used in place of the eddy diffusivity from the main TCM at the boundary 
layer top to further calculate buoyancy flux used to obtain the time rate of change in TKE 
and for turbulent mixing of scalar variables. 

Special considerations were given to buoyancy flux just below the boundary layer 
top for coarse resolution mesoscale models. It is well recognized that cloud-top radiative 
cooling is a major forcing of turbulence in the cloudy boundary layers (e.g., Nicholls 
1985, Lilly 1968). This radiative cooling, normally occurring in a layer of less than 100 
m near the cloud top, cannot be well represented in general mesoscale models with a 
vertical resolution of 100 m or more. Therefore, it is advantageous to explicitly represent 
in-cloud radiative cooling in the buoyancy flux in order to capture the forcing of 
turbulence by cloud-top longwave cooling even at lower resolutions. Equation 5.6 below 
is used for buoyancy flux at the boundary layer top following Grenier and Bretherton 
( 2001 ): 


pi §_ 

b o vl 


§ M+l f 


CApny 


M + 1 +W ' 6 vl 


(5.6) 


where P,' is the buoyancy production term in the TKE budget equation, S M1 'F is the 

- ~i 

difference of the net radiative flux at the interface between cloud and clear, and w'6 v is 
the entrainment flux. 

In the single column model the effects of entrainment can be easily seen. 
Sensitivity tests were made on variables involved in the explicit entrainment 
parameterization. We first examined the sensitivity of the entrainment velocity 
formulation to the choice of the TKE value (e) in Equation 5.1. Strictly speaking, one 
should use the TKE at the diagnosed boundary layer top. However, because of the 
discontinuity involved when the boundary layer top jumps from one grid level to the 
next, using the TKE at the boundary layer top may introduce discontinuities in the results. 
Alternatively, we use the layer-averaged TKE in Equation 5.1 which reduces the 
discontinuity and reduces the tendency to produce unrealistically large entrainment 
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velocities when the inversion height is changing grid levels. Using the layer average 
reduces noise in the model and the growth of the boundary layer and the boundary layer 
thermodynamics remain about the same. In the full 3-dimensional model, it was found 
that the method for calculating the mean TKE would produce spurious values, so the 
results presented in this chapter were obtained using the TKE at one level below the 
inversion, the magnitude of which was found to be comparable to the mean TKE in the 
SCM. 

Several formulations for calculating the entrainment velocity and introducing it 
into the mixing were tested. The results presented here use the TKE at one level below 
the inversion (kj nve rsion+l) in the calculation of the entrainment velocity, which is then 
introduced into the TKE budget equation as the eddy-mixing coefficient (Kh) at the 
inversion level. Analysis of the average TKE for the boundary layer compared to the 
value at ki nve rsion+l indicates that this value is approximately equivalent to using the layer- 
integrated average TKE in Equation 5.1 and is easier to implement. The scalar mixing 
routines interpolate the Kh profile to flux levels before calculating the vertical gradient of 
the mean scalar quantity and the resulting fluxes. As this tends to reduce the effect of the 
explicit entrainment parameterization, which is now applied at the mean grid level rather 
than the flux level in the restricted inversion, the diagnosed K H jnv was applied directly at 
the flux level in the mixing routines. 

EEP implemented in COAMPS™ was tested in a SCM version of COAMPS™ 
before it was run in 3-D mode. The SCM testing allows us to understand the behavior of 
the new parameterization implementation. Figure 5.3 shows the initial/upper air 
conditions used to force the SCM simulations. Here, our SCM testing did not incorporate 
the effects of large-scale subsidence. Because large-scale subsidence is not included, an 
equilibrium solution is not expected even though the sea surface temperature remains 
constant (T sst =290 K) and boundary layer growth continues through the simulation 
period. 

The effects of vertical resolution were examined in the SCM with control and 

EEP tests done for 30 vertical levels at 100 m spacing in the boundary layer and 40 levels 

at 50 m spacing. We found that increased resolution resulted in faster boundary layer 
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growth and reduced liquid water content in both the control and EEP simulations (not 
shown), indicating that the implicit method approaches the explicit method at higher 
vertical resolutions. However very fine vertical resolution makes this sensitivity testing 
less applicable to development for use in operational mesoscale models. 

The results from SCM testing are shown in Figures 5.4 to 5.10. Figure 5.4 shows 
a comparison of the cloud water between the control run and the EEP simulation. In 
addition to the elevated cloud layer, the EEP run produces less cloud water compared to 
the control run. Entra i n m ent velocity and the associated fluxes are shown in Figure 5.5 
and the cloud top jumps in liquid potential temperature (A0i) and total water (Aq t ) are 
shown in Figure 5.6. Here mean entrainment velocities (Figure 5.5a) are on the order of 
0.5 to 3.0 cm s' 1 while larger values of 3 to 10 cm s' 1 are observed early in the simulation 
and just before the diagnosed inversion jumps to the next level. This range of entrainment 
velocities is generally larger than observed values of 0.4 cm s' 1 for similar conditions 
during the Dynamics and Chemistry of Marine Stratocumulus experiment (DYCOMS) 
(Stevens et al. 2003a) but, given the large uncertainty in field measurements of 
entrainment velocity, are still acceptable. The entrainment fluxes in the first 10 hours of 
cloud formation are large, with average values to about -75 W nf for liquid water 
potential temperature flux and 100 W m' for total water flux. Since Figure 5.6 shows 
reasonable values of the cloud top jump conditions in A0i and Aq t , this overestimate of 
entrainment flux is likely a result of strong turbulence in the early hours of cloud 
formations (Figure 5.7). The cloud top jump conditions in Figure 5.6 show a pattern of 
sharp increases in the jump strength followed by gradual decrease. This is because in a 
growing boundary layer, the boundary layer top is identified on a discrete grid level, at 
which the jump gradually decreases in magnitude due to entrainment mixing until the 
boundary layer top is diagnosed at the next higher level. 

Development of the boundary layer based on bulk Richardson number in the 
original CO AMPS and the restricted inversion height are shown in Figure 5.8. In both 
methods of defining the boundary layer height, the EEP simulation produces a higher 
boundary layer and higher cloud layer (Figure 5.4). In Figure 5.8a, after 34 hours of 
simulation we find a stable layer below the cloud base when the simulated boundary layer 
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is apparently decoupled, most evident in q t profile shown in Figure 5.8d at hour 35 and 
45. The decoupling is also represented by the double-peak feature in the TKE profile at 
hour 45 (Figure 5.9e). When entrainment is explicitly represented, the boundary layer 
tends to be better mixed in spite of a higher boundary layer. On the other hand, the 
boundary height diagnosed based on liquid water and/or 0i gradient better represents the 
boundary layer height. The decoupling is more evident in the q t gradient, although a 
slight change in the 0i gradient, indicative of a weak sub-cloud stable layer, is also 
present. 

The most significant effect of the EEP is to lift the cloud base height and increase 
the growth rate of the boundary layer height as can be seen by comparing Figures 5.9 a) 
and b) to Figures 5.10 a) and b). In addition to a faster growing boundary layer, the EEP 
simulation results in a higher and thinner cloud with liquid water confined to 
approximately one or two grid levels, while in the control run the cloud layer spreads 
over two or more grid levels. In general, the magnitude of buoyancy flux is reduced in 
the EEP run compared to the control run as seen in Figures 5.9 and 5.10 f) and g), 
however the buoyancy flux at the flux level immediately above the grid level with cloud 
water is enhanced. This is mainly due to our modification to buoyancy flux at the 
interface to explicitly anchor the long wave radiative cooling at the interface. This 
compensates for the negative heat flux introduced by enhanced entrainment so that the 
TKE in the boundary layer is maintained at approximately the same level in both 
simulations. 

Overall, the testing of the EEP in the SCM COAMPS seems to yield reasonable 
results consistent with enhanced entrainment for STCU boundary layer. Our sensitivity 
tests suggest that the growth of the boundary layer is rather sensitive to the vertical grid 
resolution, especially in our SCM testing without large scale subsidence while the 
sensitivity to choice of the velocity scale (or TKE) used in the model is less significant. 
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Figure 5.3. Vertical profiles of a) potential temperature and b) water vapor used 
to initialize the SCM CO AMPS. 
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Figure 5.4. Time evolution of the vertical profiles of cloud liquid water mixing 
ratio (g kg 1 ) from SCM tests for a) control run (CNTRL) and b) simulation with 
explicit entrainment parameterization (EEP). In both cases clouds begin to develop 
at hour 10, however the presence of a cloud layer triggers the EEP, which reduces 
the LWC for marginal cloud and delays the formation of significant cloud until 10 
hours later in the EEP simulation. 
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Figure 5.5. Entrainment velocity (cm s' 1 ) a) and entrainment fluxes (W m' 2 ) of 
liquid water potential temperature b) and total water c) from the SCM EEP 
simulation. 
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Figure 5.6. Jump condition in a) liquid potential temperature (K) and b) total 
water mixing ratio (g kg' 1 ) for the EEP SCM. 
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Figure 5.7. Time evolution of the vertically averaged square root of TKE (m s' 


*)• 


86 






































Time (hr) 

Figure 5.8. Diagnosed boundary layer height for the a) reconstructed inversion 
diagnosed using the Richardson number and b) boundary layer height restricted to 
flux level for SCM control and explicit entrainment parameterization (EEP) 
simulations. The Richardson number based boundary layer height is interpolated 
between the grid level at the inversion base and the next highest model level based 
on the strength of the thermodynamic gradient. The boundary layer height is only 
output when cloud is present and the jump condition for entrainment is met as 
discussed in the text, which accounts for gaps in the data just after the inversion 
level changes. 
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Figure 5.9. Time evolution (hours) of vertical profiles of mean and turbulence 
properties from the control run. (a) cloud liquid water content; (b) net radiative 
flux; (c) liquid water potential temperature; (d) total water specific humidity; (e) 
turbulence kinetic energy; (f) potential temperature flux; (g) liquid water potential 
temperature flux; and (h) total water flux. 
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Figure 5.10. Same as in Figure 5.9., except for EEP simulation. 
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VI. COAMPS™ SIMULATIONS WITH EXPLICIT 
ENTRAINMENT PARAMETERIZATION 


Two phenomena that were frequently observed in the inner coastal zone during 
DECS were the onset and evolution of the MLLCJ and variations in the diurnal 
dissipation, reformation, structure and inland extent of the coastal stratus. As both of 
these are closely related to the entrainment process, model simulations of these events 
with and without an explicit entrainment parameterization make an interesting case study. 
Much of the current research (e.g. Bretherton et al. 1999, Moeng et al. 1999, Moeng 
2000, Stevens et al. 2003a and 2003b) is primarily in nocturnal, buoyancy-driven marine 
boundary layers and stratocumulus transition while our study focuses on the daytime 
coastal zone structure as compared to observations. 


A. CASE STUDY CONTROL RUN SIMULATIONS 

In several of the flights during DECS a rapid rise in cloud base and lowering in 
cloud and inversion height was observed approximately 100 km from the coast. This rise 
in cloud base is consistent with decoupling and is supported by estimates of the lifting 
condensation level (LCL) based on mean parcel thermodynamic characteristics measured 
on the low level flight legs as reported in Kalogiros and Wang (2001). The LCL is the 
expected cloud base height for a parcel lifted adiabatically from the surface and is much 
lower than the observed cloud base, indicating that the cloud layer is separate from or 
decoupled from the surface fluxes. Decoupling frequently occurs when the cloud layer is 
heated due to solar radiation or increased entrainment of warm, dry inversion air or when 
the sub-cloud layer cools due to re-evaporation of drizzle, creating stable stratification. 
The rapid increase in cloud base may be indicative of decoupling, while the rapid 
lowering of cloud top and inversion height likely indicates increased subsidence. Further 
investigation from the simulation indicates this subsidence has a distinct mesoscale 
structure and is likely associated with the large-scale flow interacting with the terrain. 
The SCM simulations with EEP showed less tendency than the control run to decouple as 
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the boundary layer grew. This is likely due to increased in-cloud TKE and more turbulent 
transport between the cloud and sub-cloud layer in the EEP cases. 

A 48-hour simulation implementing the explicit entrainment parameterization was 
generated for the period from 0000 UTC July 9 to 2400 UTC July 10, 1999. This period 
was chosen for the initial testing of EEP in the full three dimensional COAMPS™ 
because the synoptic conditions showed the typical summertime transition from a fully 
cloud covered MABL with minimal LLCJ signature on July 9 to moderately strong jet 
conditions and offshore clearing on July 10. Figure 6.1 depicts the evolution of the cloud 
fields for this period. In Figure 6.1.a) at 0000 UTC the region is characterized by clearing 
to the north and northwest with solid stratus to the south and offshore. Overnight the 
winds relax and a solid persistent stratiform cloud deck develops to the coast with 
nocturnal coastal and valley fog. By 0000 UTC on July 10 (late afternoon on July 9 th 
local time) a narrow band of clearing can be seen in the satellite images in Figure 6.1.d) 
downstream of Cape Mendocino, Pigeon Point, and Point Sur indicative of the initial 
development of a coastal jet. During July 10 th broad clearing occurs offshore, especially 
downstream of Cape Mendocino, indicative of a fully developed LLCJ while the CTBL 
persists in the inner coastal between Cape Mendocino and the California Bight. The cloud 
field on July 10 th is more cumuliform or cellular in texture than on July 9 th . 

Figure 6.2 depicts the LWP evolution in the control run simulations. The control 
run produces a solid stratiform deck near the coast which matches well with the satellite 
observations but initiates the offshore clearing 6-12 hours too early. The boundary layer 
height diagnosed from the 0i gradient in Figure 6.3 shows higher boundary layer heights 
to the west and south to a maximum of 260 m with very low boundary layer heights of 
155 m or less within 50 km of the coast. During the night (1200 UTC) the boundary layer 
is seen to lower over land, although little change is seen over water except for the 
development by 24 hours of a zone of lower heights offshore between 37 and 38 N, 
downstream of Cape Mendocino. This correlates well with both the surface wind speed 
and the wind speed at the inversion depicted in Figures 6.4 and 6.5. Although the timing 
is about 12 hours too early, the clearing and increased winds in the simulation seem to 
correspond well with satellite data depicting cloud dissipation and ship and buoy 
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observations depicting an increase in surface winds from 15 to 20 knots. Despite the 
differences in the synoptic scale offshore cloud field between the observations and the 
control run, our primary focus is in the coastal zone and the mesoscale structure of the 
boundary layer. In this region the simulation does reasonably well and can be used to 
compare with the observations and the results of explicitly parameterizing entrainment. 
It is expected that the differences between the control and EEP simulations will primarily 
be in the fine scale structure of the CTBL, which will be examined in the next section. 
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Figure 6.1: Visual satellite image and buoy observations for a) 0000, b) 1300, c) 
1800, and d) 2400 UTC July 9, 1999 and e) 1400 and f) 1800 UTC July 10, 1999 
depicting solid cloud cover and light surface winds transitioning to clearing 
offshore and stronger surface winds. 
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Figure 6.2. Control run forecast LWP (kg m' 2 ) for a) 0000, b) 1200, c) 1800, and 
d) 2400 UTC July 9, 1999 from an initial analysis at 0000 UTC depicting solid 
cloud cover in the inner coastal region with clearing propagating southward west 
of 124N. 
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Figure 6.3. Control run inversion height (m) based on 0i gradient for a) 0001, b) 
1200, c) 1800, and d) 2400 UTC July 9, 1999 depicting higher boundary layer 
heights offshore and nocturnal lowering over land. The boundary layer height 
diagnostic based on virtual potential temperature gradient is not produced at 
analysis time so the 0001 UTC field is substituted. 


96 








































Control Run 

199907091 UTC Windspeed at the Surface (m s “1^99070912 UTC Windspeed at the Surface (m s' 1 )) 



1999070918 UTC Windspeed at the Surface (m sT^J)9070924 UTC Windspeed at the Surface (ms "*)) 




Figure 6.4. Control run surface wind speed (m s’ 1 ) and direction (arrows) for a) 
0000, b) 1200, c) 1800, and d) 2400 UTC July 9, 1999 depicting generally light 
winds along the coast increasing offshore. 


97 














































Control Run 


199907091 UTC Windspeed at Inversion (m s 1 )1999070912 UTC Windspeed at Inversion (m s "’j 




1999070918 UTC Windspeed at Inversion (m s ^999070924 UTC Windspeed at Inversion (m s "'j 




Figure 6.5. Control run windspeed (m s' 1 ) and direction (arrows) at the inversion 
height for a) 0000, b) 1200, c) 1800, and d) 2400 UTC July 9, 1999 depicting 
generally development of a coastal low level jet offshore. 
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B. CASE STUDY SIMULATIONS USING AN EXPLICIT ENTRAINMENT 

PARAMETERIZATION 

1. Effects Of Explicit Entrainment On Simulating The Coastal Boundary 
Layer 

Although the effect on total boundary layer growth was much smaller in the 3-D 
model than in the SCM, there was some modest improvement in the mean structure as 
compared to the observations. The LWP evolution in Figure 6.6 shows a similar trend, 
horizontal structure, and cloud edge for the EEP as in the control run, although absolute 
amounts are about 25% or 0.05 to 0.15 kg m‘ lower except very close to the coast where 
both simulations show values in the 0.4 to 0.5 range. The wind speed and direction in 
Figure 6.7 and Figure 6.8 for the EEP simulations showed no noticeable difference from 
the control run, which is not unexpected as the current EEP formulation may not 
completely represent momentum mixing. 

Figure 6.9 depicts the time evolution of the boundary layer in an east-west cross- 
section from Monterey Bay, approximately corresponding to the aircraft flight track. 
Although the EEP cloud top is not noticeably higher, the LWC is reduced and the cloud 
base is elevated, which better matches the observations. The obvious exception to this is 
within 25 km of the coast, where both simulations saturate the boundary layer to the 
surface and produce fog. This is the region of highest complexity and largest variance in 
the differences between and observed inversion heights and warrants further study. 

Figure 6.10 depicts the difference in the height of the inversion base based on the 
liquid potential temperature gradient for the control and EEP simulations. In general the 
boundary layer is 50-100 m higher in the EEP simulations, especially between 122 and 
124 degrees west, with a less noticeable effect very close to the coast where the boundary 
layer is lowest and SSTs are coldest. The exception to this is very close to the coast in 
the region of the MLLCJ, where the EEP actually produces lower boundary layer heights. 
Indications are that the current EEP underestimates entrainment in shear-dominant 
conditions. The role of shear in entrainment efficiency across the boundary layer top has 
been minimally explored in the literature and especially warrants further study. 
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With Explicit Entrainment 

199907090 UTC Liquid Water Path (kg m -2 ) 1999070912 UTC Liquid Water Path (kg m -2 ) 




1999070918 UTC Liquid Water Path (kg m 2 ) 1999070924 UTC Liquid Water Path (kg m“ 2 ) 




Figure 6.6. EEP simulation forecast LWP (kg m' 2 ) for a) 0000, b) 1200, c) 1800, 
and d) 2400 UTC July 9, 1999 depicting solid cloud cover in the inner coastal 
region with clearing propagating southward west of 124 N. Note cloud edge is 
similar to the control run but the LWP values are generally lower. 
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With Explicit Entrainment 

199907091 UTC Windspeed at the Surface (m s ”1^99070912 UTC Windspeed at the Surface (m s” 1 ) 
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Figure 6.7. EEP simulation surface wind speed (m s' 1 ) and direction (arrows) for 
a) 0000, b) 1200, c) 1800, and d) 2400 UTC July 9, 1999 depicting generally light 
winds along the coast increasing offshore. 
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With Explicit Entrainment 
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Figure 6.8. EEP simulation windspeed (m s' 1 ) and direction (arrows) at the 
inversion height for a) 0000, b) 1200, c) 1800, and d) 2400 UTC July 9, 1999 
depicting the development of a coastal low level jet offshore. 
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Figure 6.9: Time series contours of LWC (g kg' 1 ) (shaded) and 0i (solid lines) for 
the control and EEP simulation along an east-west cross-section from Monterey 
Bay. Note the reduced LWC and higher cloud base in the EEP simulations. 
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EEP-CNTRL 


199907091 UTC Inversion Height (m) 1999070912 UTC Inversion Height (m) 
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Figure 6.10. Difference between the diagnosed inversion height (m) in the control 
run and EEP simulation for a) 0001, b) 1200, c) 1800, and d) 2400 UTC July 9, 
1999 from an initialization at 0000 UTC. Note the lower inversion height near the 
coast for the EEP during the early part of the simulation and the generally higher 
heights later on in the forecast. 
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2. Changes in the Resulting Turbulent Flux Profiles and Boundary Layer 
Structure 

In order to better analyze the detailed effects of EEP as compared to the idealized 
SCM results, individual profiles were generated for the control and EEP simulations 
corresponding to the observed aircraft soundings. Figure 6.11 and 6.12 are two 
representative profiles of the mean boundary layer structure and the turbulent fluxes. The 
general wind and thermodynamic structure of the boundary layer was not seen to change 
nearly as much in the three dimensional model with EEP as in the SCM, possibly due to 
the effects of horizontal advection between cloudy and clear grid points. This is not 
unexpected as changes to the turbulence parameterization generally have the greatest 
effect in the fine scale and on long term equilibrium or climate scale. The mean boundary 
layer structure and mesoscale-to-synoptic variability respond much more strongly to the 
synoptic forcing and changes to the turbulence and in-cloud turbulent fluxes in a strongly 
subsiding environment are damped. In the three dimensional model, the EEP simulation 
depicted the same general trends and direction as in the SCM tests in that the cloud liquid 
water was reduced, the inversion was modestly higher by 1-2 grid levels (-100 m), and 
the TKE and turbulent fluxes were greater at cloud top and more negative just above 
cloud top. A sharpening of the potential temperature gradient at cloud top was also 
observed which was less noticeable in the SCM. As discussed in Paluch and Lenschow 
(1991), this sharpening of the inversion structure is characteristic of the STBL as a result 
of in-cloud TKE primarily generated by cloud-top radiative cooling responding to mean 
subsidence above the cloud. The SCM tests did not account for the presence of large 
scale subsidence acting against the growth of the turbulent boundary layer. 

The TKE near the surface was similar or slightly reduced in the EEP cases but was larger 
at cloud top. The buoyancy flux was also larger in the cloud layer and larger in the 
negative sense above the cloud which is consistent with the SCM tests. 
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Figure 6.11. Model vertical profiles (diamonds) for control run (solid line) and 
EEP (dashed line) and observed sounding leg (dark thick line) for thermodynamic 
and turbulence structure in the coastal boundary layer. 
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Figure 6.12. As in Figure 6.11 except for another location. 
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VII. SUMMARY, CONCLUSIONS, AND RECOMMENDATIONS 


A. SUMMARY AND CONCLUSIONS 

This thesis research focuses on improving the understanding and representation of 
physical processes in the marine stratocumulus-topped boundary layer and their 
interaction with the coastal flow field. Careful analysis of the observed summertime 
cloud-topped boundary layer structure from the DECS field study supported our general 
conceptual model of the coastal zone including such features as the increase in boundary 
layer height southward and away from the coast and a sharp discontinuity in the 
temperature and moisture fields at cloud top. A general synoptic pattern was observed of 
uniform fog and stratus with a diurnal signature over land to within a few kilometers of 
the coast, followed by increasing low level winds and clearing over a 2-3 day period. It 
was observed that increases in the coast-parallel wind speed at the boundary layer top 
were coincident with a distinct mesoscale variability in the horizontal structure and 
diurnal dissipation of the cloud, especially within 100 km of the coast. This has been 
attributed to flow interaction with the terrain creating velocity minima, compression 
bulges, higher inversion heights and more persistent cloudiness upstream of points and 
expansion fans, clearing, and accelerations downstream of points. In general, the 
observations support this model, however it was noted that in the cloud topped boundary 
layer, the velocity maximum was frequently below the inversion, or the boundary layer 
velocities were nearly uniform with height to a minimum above the inversion. This is 
somewhat different than the expected structure of the low level coastal jet and warrants 
further study. In the clear boundary layer, the structure was more classic, with a well 
defined maximum in velocity at or just above the inversion height. 

In the control run simulations COAMPS™ was used to simulate the atmospheric 
structure for the entire DECS period in order to better understand the observed patterns 
and assess the forecast model’s ability to represent them. The model winds, cloud fields 
and boundary layer structure matched the observations quite well in the synoptic and 
diurnal evolution and horizontal variability, with the exception of some consistent biases 
in the inversion height and cloud liquid water content. More detailed analyses revealed 
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that the cloud-topped marine boundary layer height in COAMPS™ is, on average, under¬ 
estimated by about 200 m and the liquid water content is over-estimated by about 25%, or 
about 0.15 g kg' 1 . Furthermore, the model matched the observations in the clear marine 
boundary layer much more closely, indicating that the bias in the control run may be 
related to particular characteristics of the CTBL. These biases were not found to correlate 
well with distance from the coast, synoptic variability, or time of day, although the 
difference between the observations and the model was much higher with more spread 
close to the coast. This is expected as flow interaction with the terrain and land-sea 
interactions make this region particularly complex. The underestimate of the cloud top 
and base has greater operational impact near the coast because the boundary layer is in 
general rather shallow, resulting in COAMPS™ producing fog rather than the observed 
low stratus in the inner coastal zone. 

The correlation with cloudy conditions implies that the bias is likely due to 
misrepresentation of cloud-related processes. The high LWC and low inversion top is 
consistent with under-representation of entrainment fluxes at cloud top. Therefore we 
implemented an explicit entrainment parameterization scheme to test this hypothesis. 
Following the work of Grenier and Bretherton (2001) we developed a framework within 
COAMPS™ to diagnose the entrainment velocity across the inversion, calculate the 
entrainment fluxes and represent the entrainment effects explicitly in the forecast model. 
We tested this framework in the convective stratocumulus-topped marine boundary layer 
in the single column model and found that, in the absence of subsidence, the explicit 
entrainment fluxes caused the boundary layer to grow faster than in the control run, the 
boundary layer was better mixed, the potential temperature was higher, the cloud layer 
was thinner and the liquid water and total water content was lower. These trends were 
consistent with explicit entrainment being more efficient than the implicit method used in 
the parent turbulence closure model at entraining warm, dry air into the cloud layer from 
above the inversion. 

Some technical issues were examined in the SCM formulation and testing. These 
included several approaches for specifying the inversion height in a robust and consistent 
way. This is especially critical when the boundary layer is evolving rapidly and the 
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inversion strength and height is changing between discrete grid levels because 
inconsistent treatment results in incorrect entrainment flux and model noise. The most 
robust method was found to be direct diagnosis of the cloud top, which is simple to 
implement but restricts further study of entrainment with this EEP framework to the 
CTBL. Another critical area is the specification of the appropriate value for TKE in the 
explicit entrainment parameterization. TKE at cloud top, one level below cloud top and 
the layer averaged TKE were tested. Another possibility that was not explored is to use 
other estimates of the turbulent forcing such as the convective velocity scale (w*, Stull 
1988). Other technical issues that were addressed include the introduction of the eddy 
diffusivity mixing length consistently into the TKE prognostic equation and the turbulent 
mixing subroutines, specification of the thermodynamic structure of the sub-grid scale 
ambiguous layer that the inversion lies within, and accounting for key sub-grid scale 
processes such as evaporative enhancement of the TKE generation by mixing between 
cloudy and clear parcels and radiative flux divergence across a discontinuous boundary. 

Once a reasonably stable and consistent scheme was developed for the SCM, a 
full three dimensional case study was run to explore the impact of the EEP with realistic 
forcing and horizontal variability. The explicit parameterization was found to produce 
trends qualitatively similar to the SCM in that the inversion layer was cooled and 
moistened, the inversion sharpened, and the cloud liquid water was reduced. The surface 
buoyancy and total water fluxes were reduced somewhat but the cloud-top TKE was 
increased. In some cases, the in-cloud total water and buoyancy flux was seen to increase 
and the negative buoyancy flux above cloud top became stronger. This resulted in a 
profile that better matched the observations. The actual cloud top was highly variable 
and increased moderately in some cases and decreased in others. Overall, the net effect 
was to lift the inversion by about one grid level, or 50-100 meters and reduce the LWC so 
that the surface layer was no longer saturated. The observations indicated that the cloud 
deck was frequently 100-200 meters thick, which is difficult to accurately represent when 
the model vertical grid resolution is the same magnitude. 
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B. 


RECOMMENDATIONS FOR FURTHER STUDY 


This thesis work has established a basic framework for explicit entrainment 
implemented into COAMPS™. Based on this framework, various formulations can be 
evaluated and observational data can be used to evaluate the appropriateness of the 
associated free parameters and constants originally derived from other field studies or 
LES data. As discussed in Chapter II, several researchers have developed schemes for 
calculating the entrainment velocity that are more sophisticated than the basic local 
closure implemented here that could easily be tested in the current framework. While 
most of the current entrainment research focuses on buoyancy effects from radiative and 
evaporative forcing, a complete parameterization should also include shear effects, for 
which its expected that the MLLCJ interaction with the CTBL top and the frequently 
observed indicative cloud signature will make an excellent case study. The DECS dataset 
includes cloud radar, rawinsonde, ceilometer, RASS and profiler data as well as surface 
mesonet and aircraft data, which should facilitate the study of CTBL diurnal evolution as 
the high data rate aircraft sounding data available in many field studies can be augmented 
with excellent temporal coverage in the vertical at a fixed point. Furthermore, many of 
the efforts such as Grenier and Bretherton (2001) and Lenderink and Holtslag (2000) that 
combine explicit entrainment with a TCM were tested using single column models. As 
reported in Bretherton et al. (2004) and McCaa and Bretherton (2004), other researchers 
have examined the effects of explicit entrainment and shallow convective 
parameterizations on the climate scale mean structure of the summertime Eastern Pacific 
marine CTBL, however the mesoscale variability and vertical structure in the inner 
coastal is still relatively under studied. Implementing this approach in a three dimensional 
model with full physics will allow us to better understand the mesoscale effects of the 
three dimensional flow in the inner coastal zone. The resulting parameterization should 
be tested across a broad spectrum of regimes besides the summertime coastal marine 
CTBL before implementation in operational CO AMPS™. 
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